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Preface. 



Scientific Aspects. 

The "real space renormalization group" concept unites many of the reasons for which I like 
mathematical physics. It is a key tool (and shall remain so) in the explanation both of fundamental 
laws of Nature and human scale physics. The basic idea is accessible to high school students and, 
at the same time, deep enough to put forward some of the hardest conceptual problems we are 
facing. 

The bet for the real space renormalization group is at the heart of a certain world view as 
something fundamentally intelligible. This new physics aspires to be a constructions game: blocks, 
cells, graphs... discrete geometry and analysis as the basis of the explanation of magnetism, 
elementary particle physics, chemistry, cosmology or turbulence. This view would result too naive 
if the scales interaction and fractal geometry were not taken into account. They endow this view 
with the necessary complexity. 

Computers are the new tool, and programming the new algebra. But intelligence must stay as 
the guide, and we should not allow the tools to become the crutches of thought. I.e.: the ideas must 
be clear, every step in a computation must be deeply understood and meaningful ("never make 
a calculation whose result you do not know"). When we are able to write an efficient algorithm 
to predict a physical phenomenon, the reason is that we understand it well. When a problem is 
solved through computational "brute force" , physics have failed. 

Methodological aspects. 

All computations in this thesis have been carried out with our own programs or with free 
software. The usage of propietary programs is, for the present author, incompatible with scientific 
work. It is equivalent to the usage of mathematical formulae whose proof is not accessible to us. 
Programs have been written under LINUX in the C and C++ languages (with Emacs) for the GNU 
compilers. Libraries were written (see appendix B) for the usage of matrices, graphs, screen (Xll) 
and printer (PostScript) graphics managers... which shall be released to the public domain in 
due time. This thesis is written in pure TgX, and its format was designed by the author. 

Internet is changing the social aspects of physics. It has already removed many spatial barriers 
to scientific communication, and is now weakening the most powerful institutional burden: scientific 
journals. As a tiny contribution, all our papers have been sent to the Los Alamos, including this 
thesis. 
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Synopsis. 

1.1. New techniques for new questions. 

1.2. The scale axis and the renormahzation group. 

1.3. Blocking physics. 

1.4. Brief history of the renormahzation group. 

1.5. Plan of this work. 

1.6. Bibliography. 

The present chapter consists of a non-technical introduction to the Real Space Renormahzation 
Group (RSRG) and the associated mathematical-physics concepts. There are no new results, 
although the exposition is original. The purpose is merely to make the whole work as self-contained 
and accessible to the non-specialist as possible. 



1.1. New Techniques for New Questions. 

It may be argued that the biggest renewal movement in theoretical physics since the middle of 
the past century has been the relevance of the analysis of complex systems whose basic rules are 
known. Nowadays, almost all resarch fields involve the study of systems with many (maybe infinite) 
degrees of freedom. A fast survey of the key problems range from the large-scale structure of the 
Universe [KT go], [PEB 8o], high critical temperature superconductivity [GMSV 95], hydrodynamical 
turbulence [FR95], [MY 65], surface growing phenomena [BS95] or quark-gluon plasma physics 
[DGMP 97]. 

It is probably excessive to assert that complex phenomena physics involves a new scientific 
Weltanschauung^. But, nevertheless, it is important to acknowledge the truth contained in the 
proposition: 

• There are new methods: doing physics involves different jobs than before. Maybe the most 
From German: world view. 
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important one is the extensive usage of computers. 

• There is a new series of "paradigmatic examples" about what is an appropriate answer to 
a problem. We may cite the solution to the 2D Ising model by L. Onsager [FEY 72], the 
renormalization group analysis of the Kondo problem by K.G. Wilson [WIL 75], the turbulence 
theory of A.N. Kolmogorov [K 41] or the application of the density matrix renormalization 
group to the problem of the particle in a box [WHI 98]. 

• The questions considered to be relevant have also changed. Thus, e.g., the precise prediction of 
the power law according to which the specific heat of a ferromagnet diverges as it approaches 
Curie's point is considered to be a great success even if the position of that temperature is not 
quite accurate [BDFN93]. 

The massive usage of computers has had a great influence. A problem is considered to be well 
solved when a fast and correct computational algorithm has been provided which allows the precise 
calculation of physical observables. A greater importance is being assigned to discrete methods in 
physics, which were almost completely forlorn. Although traditional mathematical analysis was 
developed in parallel with discrete analysis-^, it was abandoned during the XVIII and XIX centuries 
but for a few geniuses such as Leonhard Euler or George Boole [GKP 89]. 

Some of the most important physicists of the era previous to this stage, to which they greatly 
contributed (such as, e.g. L.D. Landau), were against the expansion in the usage of computers. 
In a certain sense, they were right. During the technological explosion which, since the 70's, took 
computers to the physics departments all over the world, the physical and mathematical insight 
was abandoned to some extent for the mere computer-aided exploration of solutions to problems. 
The incredible boost in the capacity of the computers allowed many researchers to neglect their 
algorithms. Many problems were out of reach of their computational power, but the number of 
those which might be tackled was so big that this did not constitute a serious problem. 

But a parallel movement was being started. The idea that the search for computer algorithms 
to solve problems in physics was as respectable as it had been the search for solution methods for 
differential equations, was getting mature. 



1.2. The Scale-Axis and the Renormalization Group. 



The Scale-Axis. 

The new problems in physics are characterized by the great number of degrees of freedom. 
Their objects of study are extended systems, whether spatial and/or temporally, which have a 
basic scale of actuation at which its dynamical rules are simply determined ^. 

Unfortunately (or maybe not so) the observation scale is always much greater than the basic 
scale, and an amount of non-trivial phenomena appear which are believed to be direct consequence 
of the physics at the basic scale. The search for these basic rules belongs to another "scientific 



Newton's interpolation formula is considered the discrete analogue of Taylor's expansion, and both are not 
too distant in time [KLI 72]. 

The contents of this section can be checked in any of the basic books on RSRG [BDFN 93], [GOL 93], [CFP 
92]. Nevertheless, the presentation is new and the terms basic scale and observational frame are original. 
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program" , which was preponderant in the first half of the past century. 

No observational device can have infinite resolution or infinite range. Therefore, each descrip- 
tion of a physical system must always be labelled by a minimum scale (the grain-size of the film) 
and a maximum scale (the plate-size). In the quantum field theory nomenclature, these parameters 
are respectively called "ultraviolet-cutoff" (UV) and "infrared-cutoff" (IR). 

Let us, therefore, define an observational frame to be a reference frame (whether Lorentz or 
Galilei), along with the specification of the ranged scales. The analogy with a magnifying glass 
may be interesting [GAW 96]. 



Galaxy clusters 
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Bacteria 



IR-cutoff 
Obs. frame > 
UV-cutoff 



Atoms 



Quarks 



Planck scale 



Figure 1. An observational frame which is carefully placed for the examination of atomic nuclei. 



The fundamental frame shall be the one which contains the scales immediately above and below 
the basic scale. A translation in the scale-axis is, therefore, a movement of the observational frame 

as a whole. This operation is usTially termed a Renormahzation Group Transformation (RGT). 

When the basic scale has long been surpassed, no new physics is supposed to appear. Thus, 
the dynamics is simply "transported" along the scale-axis from its focus by a set of propagators. 



The Renormalization Group. 

We might characterize the dynamics seen at an observational frame difTcrcnt from the funda- 
mental one as an "effective" or ^^renormalized" dynamics. This new dynamics might have the same 
form as the original one, and diflfer only in the values of some parameters. In general, it is illicit 
to make this supposition. 

Even so, we shall describe one of those textbook examples in which everything works properly: 
the Ising model [LEN20], [ISI25] which explains ferromagnetism in uniaxial magnets. It only 
contains classical spins which are either parallel or anti-parallel to a given fixed axis. It is a 
model with only two parameters at the basic scale: the coupling strength between neighbour spins 
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(attempting to order the system) and temperature (attempting to disorder it)^. 

Any physical magnitude defined at the basic scale may be "translated" along the scale-axis 
by the RGT. During its displacement, the relevance of the given element is either diminished or 
boosted. In the first case it shall be termed an irrelevant element, while in the second case it shall 
be said to be relevant^. 

Physics observed with the naked eye is so far away from the basic scale that the RGT are 
usually allowed to reach a fixed point. I.e.: advance until the relative importance of the different 
elements of the dynamics is invariant under the transformation^. 

The fixed point may yield a trivial dynamics if only one of the elements is relevant. Thus, in 
the Ising model of ferromagnetism, if only the coupling between spins is relevant, the phase shall 
be magnetized and homogeneous. On the other hand, if only temperature is relevant, then the 
phase shall be disordered. 

Critical points. 

The fixed point can be much more interesting. Being invariant under scale changes, they may 
have non-trivial fractal properties [DUP 89], [KRO 00], and in that case the system is said to be at a 
critical point. No deterministic fractals are involved (such as the Mandelbrot set or the Sierpinski 
triangle), but statistical ones. In other words: the probabilities for different configurations have 
scale invar iance properties. 

The term "fractal" was coined by Benoit Mandelbrot [MAN 82] so as to unify a big number of 
explanations to physical phenomena which were characterized by a wide range of scales for which 
either dynamics or geometry was self-similar. Certainly, the study of fixed points of the RGT is a 
powerful tool in the explanation of the origin of fractal structures in Nature. 

In the Ising model example, the magnet at Curie temperature (its critical point) presents do- 
mains with either up or down magnetization, which are themselves splashed with smaller domains 
which themselves... and so on ad scalam basicam. The system has a hyerarchy of domains, which 
makes it highly susceptible. Being all scales strongly coupled, local changes arc greatly magnified. 
Thus, a small external field may impose a global magnetization on the system quite easily. 

Critical phenomena, i.e., phenomena with a certain scale- invariance, are ubiquitous in the real 
world. Although the discussed example has been a phase transition in an equilibrium statistical 
model, there are many models far from equilibrium which present self-organized crisis^. Beyond 
continuous phase transitions (without latent heat), the self-gravitating media [PEB 80], turbulence 
[FRI95], granular media movements [BTW88], polymers in solution [GEN 72], surface roughening 
[BS95], etc. 

Moreover, it should not be surprising that the mathematics of scale invariance is an essential 
tool for the theories on the fundamental constituents of matter [COL 85], [KOG 79]. In fact, as it is 

^ As a matter of fact, there is only one parameter, which is the ratio of these two. Nevertheless, it is preferrable 
to state it in these terms. 

^ Of course, marginal elements which neither increase nor decrease may also exist, and are usually rather 
difficult to analyze. 

^ This really means to take the thermodynamic limit. 

^ The word "criticality" , albeit widely used, is not correct. The term crisis is more appropriate. 
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explained in 1.4. II., it was in that field where the term "Renormalization Group" was born. 
Critical Exponents and Universality. 

At the critical point there are many physical magnitudes which diverge. For a phase transition 
those might be e.g. susceptibility, specific heat, etc. These are physically meaningful divergences, 
and in many cases are related to the inability of euclidean measures to capture a fractal phenomenon 
(see [DUP 89] and [KRO 00]). In other terms, they diverge for the same reason as the length of a 
coast does, and in the a similar way: because a ID measure is not appropriate. The length of a 
coast only diverges when there is either no highest or no lowest scale, and similarly the magnitudes 
at a critical point only diverge for an infinite system. 

In order to characterize the magnitude of these divergences, power laws are appropriate: F oc 
t""*, where F is the observable, t is a variable which is zero at the phase transition (e.g.: reduced 
temperature) and a is what we shall call a critical exponent. 

Critical exponents may be obtained from the RGT. The procedure starts with the obtention of 
the fixed points of the transformation for a given physical system. Afterwards, the neighborhood 
of the fixed points is carefully studied. Obviously, only the relevant elements of the dynamics shall 
be important in the computation. 

Universality is almost trivial after the previous remarks. If two physical systems share all the 
relevant variables and differ only in the irrelevant ones, they must share critical exponents and 
they will enter the same universality class. Thus, in the ferromagnet example, the crystallization 
type of the lattice, the chemical aspects of the atoms, etc. happen to be irrelevant. Therefore, all 
ferromagnets have the same critical exponents. Moreover, models as different beforehand as the 
"mixing-non-mixing" phase transition between CH3OH and CeHi4, the "superfiuid-normal fiuid" 
transition for ^He and many others share the essentials of the ferromagnet transition and fall into 
the same universality class. 

Experience shows that the "essentials" for a phase transition (i.e.: the variables labelling 
universality classes) are merely the dimensionality of the configurational space and that of the 
order parameter (i.e.: the field which is zero at the disordered phase and non-zero at the ordered 
one). 



1.3. Blocking Physics. 



Spin Blocks. 

The basic idea of Kadanoff [KAD 66] is the image which is kept in mind by any practicioner of 
RG techniques: the spin blocks. Figure 2 shows a square lattice of spins with three levels of block 
composition on it. 

The fundamental operation is the composition of spins to form blocks, and the fundamental 
idea is that "a spin block behaves, to some extent, in the same way as a single spin does"^. A 
description of the system in terms of blocks instead of spins only requires an appropriate change 
in the coupling constants. 



Hence the quotation for this work. 
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Figure 2. A square lattice of 8 x 8 spins with three levels of composition of blocks superimposed. Each level 
corresponds to a RG step in Kadanoff 's approach. 

Unfortunately (or perhaps not so) that bucoHc image is impossible to implement in practice. 
In many cases the effective description of a system in terms of spin blocks is possible, but at 
the expense of introducing new interactions between them which were not present beforehand. 
These new interactions, at the next step, might require even more interactions and so on ad 
nauseam. Some theories are exactly renormalizable (among them all the interactions considered 
to be fundamental except gravitation), and many of them are approximately renormalizable. I.e.: 
they allow the application of the RG to a good precision with the inclusion of a small number of 
extra interactions. 

Once the set of interactions which is approximately closed under RG has been delimited, it is 
only needed to find which is the appropriate change in the coupling constants. A parameter space 
is established with them all, and the RG has a "realization" on it as a Bow. This flow has probably 
some fixed points which, as it was discussed in the previous section, determine the macroscopic 
physics of the system. 

In this text we shall not enter the rather interesting theory of fixed points of the RGT and 
the determination of the critical exponents from the study of the flow around them. This work is 
heading in other direction. 

How TO Make up a Block. 

In spite of the conclusions which might be drawn from the previous section, the formation of 
a block from its constituent spins is not a straight-forward procedure. In physical terms, we must 
determine which are the variables which will characterize the block. 

Let us consider an Ising-like model with two states per spin (+1 and —1) and the formation 
of blocks, as in figure 2, through the composition of four neighboring spins making up a square. If 
the block spin must be an Ising spin again, it must take a value +1 or —1. How should that value 
be chosen? 

• Through decimation: the value of the upper- left spin is chosen (for example). It seems to be 
an arbitrary method, but there are suitable techniques for its realization. 

• Through majority: the sign of the sum of all the spins is chosen. It has one inconvenient: it 
does not provide a rule when the block magnetization is null. In that case, adopt a random 
value. 

The inconvenient of the majority rule is solved if the restriction of Ising spins is lifted, and the 
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adoption of arbitrary values for the spin is allowed. In that case, the rule shall also be simple: take 
the average value of the spins of the block (or the total spin, which is equivalent modulo a scaling 
factor). 

A different way to focus the problem may provide some insight. From the four original degrees 
of freedom (see again fig. 2), RG must choose a single one. According to the majority rule, it shall 
be the mode which is homogeneous on all four spins. Therefore, the system is projected on the 
degrees of freedom which are homogeneous for scales below the block size (which becomes a new 
UV cutoff). 

Is it necessary to choose precisely that degree of freedom for the block? By no means. An 
interesting option, which shall appear to have outmost importance for the rest of this work, is the 
adoption of different degrees of freedom for each block, depending on their circumstances (e.g.: 
distance to the border, values of some inhomogeneous parameters...). Moreover, each block may 
be represented by more than one degree of freedom, allowing a more refined analysis. 



1.4. Brief History of the Real Space Renormalization Group. 



I. Prehistory 

The idea of scale invariance is old and venerable in geometry. Scaling arguments were com- 
monplace among the pithagorean school/sect (the golden section and its realizations), culminating 
with Euclid's proofs (see, e.g. [KLI72]). These ideas, held alive through the Middle Ages in a 
quasi- mystical disguise, had a great influence on the re-foundation of physics in the XVII. For ex- 
ample, Galileo [GAL 1638] makes a beautiful discussion on the scaling properties of animal sizes^. 
Notwithstanding, the mathematical tools of the rationalist era, early algebra and analysis, were 
not well suited to deal with problems with many scales. The appropriate tools may have appeared 
only recently. 

Until the end of the XIX century, scaling arguments were not common in physics. Fluid 
mechanics ([MY 65], [FRI95]) is a classical source of successful scaling hypothesis. The idea of 
Osborne Reynolds about considering velocity fluctuations in a flow as the generators of an "effective 
viscosity" is probably one of the first examples of a parameter renormalization in the modern sense. 
Moreover, the qualitative image of Lewis F. Richardson (1922) about turbulence-^" was converted 
by Andrei N. Kolmogorov (1941) [K 41] into a quantitative theory of a suprising success. 

By the beginning of the past century, acoustics and biophysics had already generated a great 
amount of such kind of arguments, although most of them were merely empirical [SCH 91]. 



E.g., if an animal doubles it size, its weight is multiplied by eight, while the section of its bones is only 
multiplied by four, rendering it unstable. 

Lewis F. Richardson was a peculiar meteorologist with a strong intuition and unorthodox expressive means. 
He paraphrased a poem by Jonathan Swift [MY 65]: 

Bigwhorls have little whorls, 

which feed on their velocity; 

and little whorls have lesser whorls, 

and so on to viscosity (in the molecular sense). 
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The last decades of the XIX century saw the rebirth of the attempts to explain human scale 
phenomena on the basis of the movement and interactions of atoms^^: kinetic theory and statistical 
physics. On the other hand, the predominant scientific philosophy was positivism, which proclaimed 
that deep understanding was not a main objective of physics, and discarded atoms as unnecessary^^. 
Albert Einstein was partly guided by scaling arguments in his search for a diffusion theory [EIN 
05] , whose importance is similar to his other two 1905 papers, since it destroyed the opposition to 
the atomic hypothesis and statistical physics. 

Statistical physics opened up a new mathematical problem: the study of many similar ele- 
ments (molecules, spins, electrons...) which create a coercive field on their neighbours. The most 
immediate ancestor of the renormalization group is probably Mean Field Theory [BDFN 93], [GOL 
93] , which is an attempt of general applicability to understand such systems by assuming the same 
behaviour for all the items and searching a self-consistent solution. In other words, each element 
must act in such a way that its behaviour gets explained by the interaction with neighboring 
elements behaving as it does 

Some of the prototypical examples (of great success) are Van der Waals theory of the liquid- 
vapour transition, Weiss theory of ferromagnetism and Hartrce-Fock theory of atomic spectra [HAR 
57] (which was extended to interacting electrons in solids [MAT 67]). 

Mean field theory predicted a phase transition for the Ising model of ferromagnetism in any 
dimension. Since the theory presupposes the absence of dissent among the spins (there are no 
fluctuations), the higher the number of neighbours of each item is, the better it works. Specifically, 
in dimension > 4 mean field yields exact results for many theories. But for unidimensional problems 
it is unable to predict that, due to geometrical scaling reasons^^, there are no phase transition. 

Mean field theory reached its most highly sophisticated expression with the ideas of L.D. Lan- 
dau, who employed it to analyze the first great statistical theory of continuous fields which served 
to model a big number of phase transitions with a mere change of its parameters. V.L. Ginzburg 
proceeded to analyze the influence of the fluctuations and to give an applicability criterion for mean 
field on Landau theory. But the mathematical theory was so intrincate that it was impossible to 
improve the approximations for the critical exponents given by mean field until the introduction 
of the RG. 

Since the end of the XIX century a series of concepts which were later known as fractal geometry 
were getting mature [MAN 92], [EDG 90]. The monster collection of classical analysis (Peano curves, 
Weierstrafi functions...) along with the insightful comments of Henri Poincare on the nature of 
non-integrability, the measure theory of Felix Hausdorff (which includes non-integer dimensions) 
joined the powerful intuition of Benoit Mandelbrot to create a geometric frame which has proved 
to be better suited than the old euclidean scheme for the understanding of phenomena with many 
scales^^. 

Wo have said rcbirtli because yet Plato's "Timaeus" was tlie first textbook on that subject [PLA 360 B.C.]. 

12 

It is said that Ludwig Boltzmann commited suicide partly because of the reject by Ernst Mach and his 
school. 

Act in such a way that your behaviour may serve as a general rule: the kantian cathegorical imperative 
provides an interesting analogy within ethics [KAN 1788]. 

"'"^ The frontier of a unidimensional region has a magnitude which is independent of its size. Thus, if a fluctuation 
is created with downwards spins in the middle of a sea of upwards spins, the cost is the same whether it is small 
or big. These fluctuations destroy any possibility of ferromagnetism (long range order) at any flnite temperature. 

Perhaps it is worth to notice Mandelbrot's comment [MAN 82] on the fact that it was UV catastrophes which 
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II. The Archaic Age. 

RG made its appearance in a rather different disguise as it is known today. An article of 
E.C.G. Stueckelberg and A. Peterman in 1953 and another one by M. Gell-Mann and F.E. Low in 
1954 [SP 53], [GML 54] opened the field, reaching maturity with the text of N.N. Bogohubov and 
D.V. Shirkov in 1959 [BOG 59]. 

These original works were dedicated to the elimination of infinites in quantum field theories 
through the "renormalization" (change, adaptation) of the physical parameters of the system. 
They were posteriorly accused of lack of insight. In words of Kenneth Wilson [WIL 75], "the worse 
feature of the standard renormalization procedure is that it is a purely mathematical technique for 
substracting out the divergent parts of integrals in the continuum limit. It gives no insight into 
the physics of the statistical continuum limit" . 

The interpretation of the RG for quantum field theories in terms of scale invariance (or its lack 
thereof) was carried out by C.G. Callan [GAL 70], K. Symanzik [SYM 70] and K.G. Wilson himself, 
around 1970, opening the path to the new era. 

The archaic era provided us with a great amount of notation, beginning with the very term 
"renormalization group" . Despite the general consensus on the inappropriateness of the term, it is 
sure that there is no general agreement about the necessity of a change or the direction it should 
take. 

III. The Golden Age. 

The main actors appear into stage around the middle sixties. In 1966, in a short-lived journal 
called Physics, Leo P. Kadanoff [KAD 66] proposed the transformation of spins into block spins for 
Ising-like models and proved that this transformation explains some (so far) empirical relations 
among the scaling exponents. 

The honour of converting the idea into an efficient computation method corresponds to Ken- 
neth G. Wilson. In 1971 he published two consecutive papers [WIL 71 A] and [WIL 71B]. The first 
one recasts Kadanoff's theory in differential form, rendering it more suitable for mathematical 
analysis. The second one anticipated wavelets [LEM 89] in more than a decade, and used them 
to analyze Ginzbiirg-Landau's model. It constituted an intermediate step between RSRG and the 
more developed Momentum Space RG. 

The early RSRG techniques were successfully applied to the 2D Ising model (see, e.g., Niemeyer 
and Van Leeuwen [NvL 73]). In 1974 Wilson solved the Kondo problem, which dealt with the effect 
of a magnetic impurity on the conduction band electrons of a metal [WIL 75]. He remarked that 
"it was the first example where the renormalization group program had been carried out in full" . 
His solution was based on the division into shells of the whole lattice around the impurity, in such 
a way that the further a shell was from the center, the nearer to the Fermi surface the electrons 
were supposed to be^^. Shells were integrated in an iterative way, starting from the center, so the 
external ones only saw the impurity sheltered by the inner shells (reminding somewhat of Gauss' 
theorem). It is important to remark that the solution required a strong amount of numerical 
computations (in 1974 terms) and it was considered to be a great success with an error in the 

killed both classical physics (fl900) and classical (smooth) mathematics (fl875). 

The usage of mixed real space and momentum space techniques is one of Wilson's great ideas. Electrons 
which are far away from the impurity only contributed to the magnetic susceptibility (his main target) if they 
were so near to the Fermi surface that their ability to be excited could compensate. 
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observables of a few percent. 

It was also in 1974 when Leo Kadanoff and his collaborators worked out the first results for 
one of the most significative problems: the 2D and 3D Ising model [KAD 75] [KH 75]. An advance 
of an extraordinary significance for the applicability of the RG to the computation of the critical 
exponents of phase transitions was the Migdal-Kadanoff transformation. It is a variational relation 
between the free energies of two lattices which differ only in the position of some bonds, as it can 
be seen in figure 3. 



Figure 3. Migdal-Kadanoff transformation. Moving the marked bonds in the given directions leave an unhnked 
site. 



The right side lattice in figure 3 has a spin site which is completely isolated, and hence its 
contribution is trivial. The free energy of the right lattice is a lower bound for the free energy of 
the primitive lattice. Although the transformation was initially developed by A.A. Migdal in 1976 
[MIG 76], it was the interpretation given by L.P. Kadanoff in terms of moving bonds which made 
it popular [KAD 76]. 

By the end of the 70's the idea of developing the RSRG as a non-perturbative computational 
method to obtain results for quantum field theories was getting mature. Some groups working in 
parallel, one in the United States ([DWY76], [DWY77], [DSW78], [DW78], [SDQW80]) and others 
in the Netherlands and Paris ([BCLW77], [PJP82]) developed similar techniques. The first one 
had as their target the low energy properties of quanlTim chromodynamics (QCD), meanwhile the 
others focused principally in condensed matter many-body theories. 

The method was baptized as the Blocks Renormalization Group (BRG). A detailed explanation 
may be found in chapter 2, but the general idea deserves to be exposed here succintly. The objective 
is to obtain variationally the lower energy states of bigger and bigger blocks using as Ansatz the 
lower energy states of the smaller ones. 

In 1982 Kenneth G. Wilson was awarded (in solitary) the Nobel prize for his RG solution 
to the Kondo problem. It is interesting to notice that by that time RSRG was going through 
an impasse, since the results of the previously stated applications to QCD and condensed matter 
physics were not as successful as it was hoped. On the other hand, big computer facilities allowed 
the implementation of massive Quantum Monte-Carlo [BH 92] algorithms to solve some of the same 
problems. Thus, RSRG was slowly abandoned. 

By the same time, the ideas of scale invariance had already invaded the theoretical analysis 
of quantum field theories (Callan-Symanzyk equations, scaling anomalies...) [COL 85], [ITZD89]. 



1.4. Brief History of the Real Space Renormalization Group. 
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Practical Monte-Carlo calculations for lattice field theories required the conceptual apparatus of 
RG to take the "statistical continuum limit" , as Wilson had foreseen. 

IV. The Scholastic Age. 

The (relative) abandonment by the groups which were mostly interested in the obtention of 
numerical results left the field ready for deeper theoretical analysis. Through the 80's a series of 
algebraic and conceptual techniques were developed which helped to understand the RSRG process. 

Touching briefly a topic which surpasses this brief introduction we shall remark the analysis of 
Belavin, Polyakov and Zamolodchikov in 1984 [BPZ 84], [ITZD 89], which issued the idea that 2D 
universality classes might be thought to be linked to representations of an algebra derived from 
RG^^ in their work on 2D Conformal Field Theories (CFT)-^^. This work triggered an avalanche 
of activity, mainly due to the interest which in a few years was unleashed by string theory^^. The 
application of CFT and, afterwards, of quantum groups to the problems of statistical physics was 
a task which, although originated at this time, was matured only through the decade of the 90's 
[GRS96]. 

V. The Industrial Age. 

Numerical methods and theoretical ideas which had "something" in common with RSRG had 
been developing in parallel with it. This proves that the main idea was mature and needed, both 
technical and theoretically. Among the most applications-oriented ideas we may cite multi-grid 
methods [PTVF 97], among the most theoretical, quantum groups or conformal field theories [GRS 
96] and among the mixed ones, wavelets theory [LEM 89] and multi-fractal analysis [FRI95]. 

Quantum Monte-Carlo (QMC) method, which had displaced RSRG by the beginning of the 
80's, had serious problems when applied to critical phenomena. Beyond the strong requirements 
on CPU time, memory and the so-called critical slow-down (which made the convergence of the 
observables despairingly slow in the vicinity of a critical point), it suffered from the sign problem 
[BE 92]. 

QMC is based on the conversion of a quantum problem into an equivalent classical equilibrium 
one through use of the Trottcr-Suzuki formula (Feynman and Hibbs [FH 65] established the analogy 
between a partition function and a path integral, while Suzuki [SUZ 76] converted it into a useful 
formula^O). The sign problem is related to this conversion in the case of fermions: the method yields 
absurd results, such as negative "ghost-like" probabilities. Although some tricks were developed, 
very scarce progress was made during the 80's. 

17 9 

Scaling exponents would be universal in the same sense as spins for different particles are, since J is a 

Casimir invariant of SU{2). The appropriate algebra is called after Virasoro. 

1 p 

Conformal transformations are all those which conserve angles locally, though not necessarily distances. Scale 
invariance along with translational, rotational and under inversion yields full conformal invariance. In 2D, con- 
formal invariance is special because, since any complex analytic function determines a conformal transformation, 
the group has infinite dimension. 

Quantum String Theory's [GSW 87] main job is to "count" all possible surfaces which start and end at given 
curves. The necessity to avoid repetition forces to take that sum modulo re-parametrizations. That led the 
practicioners to the study of CFT. 

The relation between quantum systems in dimension d and classical equilibrium systems in d -|- 1 seems to 
require a more fundamental explanation. 
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These disappointing results made some numerically inclined researchers turn their eyes again 
to RSRG. 

The success of Wilson's computation had never been repeated. Kondo's problem (and, in 
general, impurity problems) is, in some sense, special: it consists of a single source of interaction; 
each added shell is further away from it and, therefore, interacts more weakly. 

Wilson himself had proposed, in an informal talk in 1986, to try a very simple problem in 
order to understand the failure of the RSRG: the problem of a spinless particle in a ID box studied 
with quantum mechanics [WHI99]. A few people from the public, among them Steve R. White, 
undertook the project. 

This problem had the advantage of its simplicity and, moreover, it seemed to include all 
the elements which made the RSRG fail for more sophisticated problems. In effect, the blocks 
method (also known as BRG) obtained variationally the lowest energy states of the system using 
as Ansatz for the wave-functions the approximations which were obtained for the smaller blocks. 
Unfortunately (or maybe not so) the energies were incorrect by some orders of magnitude for big 
blocks, and did not scale properly. 

In 1992 S.R. White, after some unfruitful attempts^^, arrived at the satisfactory solution. The 
reason for which BRG did not work correctly with the particle in a box was that the block was 
isolated from its neighbours. White observed that the low energy states from the small blocks are, 
in many cases, quite bad bricks for building the low energy states of bigger blocks (see figure 5 in 
chapter 2). 

The best states to be chosen may not necessarily correspond to the lowest energy states. Let 
us consider, e.g., two blocks and two sites in ID as are depicted in figure 4. The ground state is 
found for the full system and it is projected on the left and right parts (left block + site, right 
block + site), providing us with new blocks (which have one more site) and which constitute good 
bricks to continue growing. 

Left Right 
Figure 4. Two blocks of several spins, and two intermediate sites which shall be swallowed by the blocks. 

The projection is carried out with the help of a density matrix, which provides the name of 
Density Matrix Renormalization Group (DMRG), profusely described at chapter 3. 

DMRG has reached numerical precision which is not only much higher than any one obtained 
with BRG, but also surpasses Wilson's solution to the Kondo problem. It has some handicaps, 
since it is based upon the left-right distinction. This makes precision diminish drastically in higher 
dimensions. 

At this point the historical exposition comes to an end in order to leave the stage for the 
explanations on our own work on this field. 



There was a successful method among the ones tried by Steve White and Reinhard Noack, but it was 
impossible to generalize. It used as "bricks" states from the blocks created with different boundary conditions. 



1.6. Bibliography. 
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1.5. Plan of this Work. 



Chapter 2 exposes a modification of tlie Blocks Renormalization Group (BRG) wliicli allows 
to analyze correctly quantum-mechanical problems both in ID and 2D. This technique is called 
Correlated Blocks Renormalization Group (CBRG). Most part of the material of this chapter is 
taken from the paper of M.A. Martm-Delgado, J. Rodriguez-Laguna and G. Sierra [MRS 95], albeit 
some ideas are exposed here for the first time. 

The third chapter describes in detail the application of the DMRG for problems in quantum 
mechanics, both for ID systems and for any kind of tree-like structure. As an application with 
technical interest, the excitonic spectrum of some macromolecules with fractal properties called 
dendrimers is obtained. Most part of the contents of this chapter is published in [MRS ooB], from 
the same authors, and the references therein. 

The standard formulation of the DMRG must necessarily distinguish between left and right, 
which is inappropriate for multidimensional systems. In order to analyze these, the Punctures 
Renormalization Group (PRG) was introduced, which is extensively described in chapter 4. The 
exposition is based on that of [MRSooA], from the former authors, but an application to the 
problem of localization in disordered excitonic systems is provided, taken from [DLMRS 01], along 
with some reflections on the possibility to extend the method to many-body problems. 

The formerly exposed methods are suitable for application to the analysis of numerical methods 
for partial differential equations. Two RSRG methods for the improvement of the efficiency of 
numerical algorithms, one based on the overlapping of cells and the other on more physical criteria, 
are developed in chapter 5. The source is the papers of A. Degenhard and J. Rodriguez-Laguna 
[DRLoiA] and [DRLoiB]. 

Each chapter has a bibliography attached. This work ends with six appendices and a glossary 
of important terms, either classical in the literature or introduced in our work. 
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Synopsis. 

• Part I. The Blocks Renormalization Group. 

2.1 The BRG for a quantum magnetism model. 
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2.8 Bibliography. 

This chapter develops in its first part the structure of technique called the Blocks Renormal- 
ization Group method (BRG) in some detail, along with the origin of its poor quantitative (and 
sometimes qualitative) results. The second part contains the original work of this chapter, which is 
the method known as Gorrelated Blocks Renormalization Group (GBRG), introduced by our group 
in 1996 [MRS 96]. It is an extension of the BRG algorithm which repairs its deficiencies and, in 
some sense, respects its spirit. 



Part I. The Blocks Renormalization Group. 

The name Blocks Renormalization Group (BRG) is associated with one of the realizations of 
the blocking idea of Leo Kadanoff [K AD 66] : the obtention of the low energy spectrum of a quantum 
system. The main advantage of the technique is that it is non-perturbative^ . Thus, states with 

It is non-perturbative on the coupling constants, but it admits a perturbative expansion on the inter-blocks 
interaction. See [GMSVgs]. 



— 17 — 



18 



2. The Correlated Blocks Renormalization Group. 



topologically non-trivial configurations and structures quite far away from the non-interacting 
states are equally available to the method. 

It was developed by some independent groups. One of them, working at SLAG, had as their 
main target the low energy properties of the theories composing the standard model, specially 
QCD [DWY76], [DWY77], [DW78], [DSW78], [SDQW80]. Others, working from the Netherlands 
and Paris, focused on models inspired in condensed matter physics, such as the Heisenberg model 
[BCLW77], [JUL 81], [PJP82]. 

The information contained in this first part of this chapter is not original to our work. It 
is included merely for completeness and to make the text as self-contained as possible. Classical 
reviews on these topics include [JUL 81], [PJP 82] or [GMSV 95]. 



2.1. BRG for a Quantum Magnetism Model. 



The Ising Model in a Transverse Field. 

The ID Ising model in a transverse field (ITF) is possibly the simplest model presenting a 
quantum phase transition [GMSV 95]. The hamiltonian is a summation of terms like this one 

H12 = -r(cT^ + cT^) - Ka| r,j>o [1] 

with a" and two of the Pauli matrices. The parameters f and J represent, respectively, an 
external magnetic field applied in the x-direction and the uniaxial coupling constant between 
neighbour spins along the z-direction. We shall assume both constants to be positive. 

The problem is non-trivial because of the competition of two effects: the magnetic field at- 
tempts to align all the spins in the x direction, meanwhile the coupling constant forces them to 
align in the z direction. The compromise is complicated because both components of the spin do 
not commute. Therefore, if the cr^-eigenstates basis is assumed, the states which are aligned along 
the X axis (due to the external field) appear to be "disordered" . 

The known relation between quantum models at dimension d and equilibrium classical models 
at dimension d -|- 1 is exemplified in the correspondence which may be established between the ID 
ITF and the 2D classical Ising model. Within it, the external magnetic field plays the role of a 
temperature (see the text by Feynman and Hibbs [FH 65] for general considerations, and the article 
by Suzuki [SUZ 76] for the explicit construction of that correspondence). 

The development of the application of the BRG to the Ising problem in a transverse field was 
carried out in [MDS 96], though we shall follow closely the exposition made at [GMSV 95]. 

The BRG applied the ITF Model. 

The main idea of the BRG procedure may be quite succintly described in a single paragraph. 
Let us consider two neighboring 1/2 spins interacting and let |-|-) and |— ) denote the eigenstates 

of a". The total system may be in four states | ), | — h), |H — ) and |++)- A hamiltonian matrix 

may be written for them and diagonalized. Let us suppose that only the low energy spectrum 
is our target. Then, the two lowest energy states may be retained and supposed to be the only 
possible states for the spins pair. To complete the analogy, we may label these states with the 
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signs "|— ')" and "|+')". Now the interaction between four real 1/2 spins may be studied as the 
interaction between two such "block spins" . 

If the system has only two sites, the eigenstates are given by: 





-) + a|++)) 


E = 


-VP + 4r2 




-+) + !+-)) 


E = 


-J 




-+)-!+-)) 


E = 


J 




— ) + 1++)) 


E = 





[2] 



where a is a parameter which depends on the coupling constant g = ]/2V, 

a(g) = [3] 

y 

Notice that the parameter a ranges from a(0) = to a(oo) = 1. The function a(g) interpolates 
between these to extreme situations: only magnetic field when a = and only spins coupling when 
a = 1. 

The first two states are now going to be assimilated as new "block spin" states I—') <— |G) and 
1+') |Ei). This operation is equivalent to a mapping from a Hilbert space of two spins (dimension 
4) into a Hilbert space of a single block spin (dimension 2). Let us define the operators: 

T : ^ 

such that 

T = |-'>(G| + |+'>(Ei| 
It = |G)(-'| + |Ei)(+'| 

The first of these operators, T, is known as truncation operator, which takes a "two spins" state and 
projects it on a "block spin" state. The second, T^, called embedding operator, takes a block spin 
state and returns a real "two spins" state. It is straightforward to observe that (since (— '|+') = 0) 

TT^ = I Identity on "block spins" space 

i.e.: the successive application of the embedding and the truncation operators yields the identity 
on the block states space. But this relation may not be inverted: 

T^T 7^ I Identity on real "two spins" space 



So, the operator T^'T, which first truncates a real "two spins" state into a block state and then 
tries to reconstruct it back, is not the identity. It is the projector on the subspace of the retained 
degrees of freedom. This shall be a recurrent idea throughout this thesis. 
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The next step is to write the hamiltonian for the system formed by two block spins. This 
hamiltonian has (hopefully) the same form as the previous one, so as the RG step may be iterated 
and reach arbitrarily long chains. 

Two blocks are placed side by side and allowed to interact via the same ITF hamiltonian. The 
"real" spins have indices 1 to 4, meanwhile the block spins have indices 1' and 2'. If j is a real 
spin index, Y shall denote the index of its associated block spin: j' = [j/2j. Figure 1 represents 
the numbering of sites and blocks graphically. 



j = 3 



Figure 1. Numbering of the sites and blocks for the BRG application on the ID-ITF. The continuous line 
corresponds to the intra-block links, and the dashed one to the inter-blocks link. 



Now let us define the block spin operators so as they fulfill relations which are equivalent to 
the ones between the old a's: 



^y\-')y = -\-')y ^y\+')y = \+')y 
^y\-')y = \+')y ^y\+')y = \-')y 



An easy algebraic manipulation leads us to 



JafV = , =af, 
' 72(1 + a2) ' 



These transformed operators shall enable us to re-write the hamiltonian for each block in its 
own terms. But the inter-blocks part (the dashed line in figure 1) requires a new operator: 

Tcf j o"2i + 1 = 2[l^-\'^'^'^) + ^ (Inter-blocks hnk) 

The most interesting conclusion of these equations is that, in fact, the hamiltonian for the 
renormalized states has the same form as the original one: 

-mN(r,j)Tt = AE + HN/2(r',r) 

with AE as a change in the "zero energy level" and two renormalized values for f and J. The three 
parameters are given by: 
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l-g^ J (l + a)^ 
(l + a)2 ^2(l + a2)J 



(l + a)2 '2(l + a2) ' '2(l + a2) 

We should remark that the exphcit form of the function a(g) was not used in the derivation of 
the renormahzation group transformation. Thus, the whole function may be used as a variational 
parameter with respect to which the energy may be minimized. The conditions which must be 
imposed to a(g) on physical grounds are its positiveness and the "boundary conditions" a(0) = 
and a(oo) = 1 (check [2] and the discussion immediately below it). We shall say that each function 
a(g) fulfilling these relations constitute an RG-prescription. 

The properties of the system depend only on the parameter g. Thus, it is appropriate to ask 
about its flow. From the previous equations it is easy to obtain 

g -'^(9)=2r-2i-a(g)2-a(l-a(g))2 ^4] 

For any function a(g) the previous transformation has three fixed points: g = 0, g = oo and a 
nontrivial intermediate one, which we shall denote by g^- The two trivial fixed points correspond 
to the ordered and disordered phases, respectively, and are both attractive. On the other hand, the 
nontrivial Qc is repulsive. The exact value of gc is 1/2, and corresponds to the same universality 
class as the 2D classical Ising model, if (g — gc) replaces the reduced temperature (T — Tc)/Tc (i.e.: 
as the variable which measures the "distance" to the critical point). 

From equation [4] it may be deduced that the value of Qc = a(gc) is a universal function. 
Whatever the function a(g), 

3 + 2gc 

Since ac € [0, 1], this equation implies that gc < 1/2, inequality which saturates at the exact 
value. 

The value of gc for a given prescription (i.e.: an specification for a(g)) can be found by seeking 
the intersection between a(g) and the function 



f(9) = { 



f 271^+29-1 

3 + 2g - y - / 

if g > 1/2 



This analysis is usually performed numerically. As an example, figure 2 shows both f(g) and 
a(g) for the prescription we have found previously. Numerically, we find gc ~ 0.39. 

Despite the numerical nature of the final part of the procedure, much can be said from merely 
qualitative considerations. Let us consider, as it is a common practice in quantum field theory, the 
function |3(g) = R(g) — g. If this function is positive, the coupling constant increases. The zeroes 
of this function, shown in figure 3, signal fixed points of the RG transformation. 

The function (3(g) fully characterizes the RG-fiow of the theory. In this case its description is 
rather straightforward, since the parameter space is ID. As it was stated above, the positivity of 
the function states the growth of g. Following that rule, the RG-flow is shown in figure 4. 
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0.35 

0,3 
0.25 

0.2 
0.15 

0.1 
0.05 


0.1 0.2 0.3 0.4 0.5 0.6 

Figure 2. RG prescription function (green line) a{g) and universal function f(g) (red one). Their intersection 
yields the value of gc for this prescription, which is approximately 0.39. 

0.4 
0,35 

0.3 
0.25 

0.2 
0.15 

0.1 
0-05 


-0,05 

-0.1 

Figure 3. Function f3{g), whose zero and infinite values correspond to fixed points. Concretely, the zero around 
g ^ 0.39 is the critical point of the theory. 

1 i — I i 1 — i ^ — i 1— ^ — ► 1 ► g 

0.3 Qc 0.5 

Figure 4, RG-flow chart of the Ising model in a transverse field according to our BRG prescription. Notice the 
critical point gc ~ 0,39, necessarily repulsive. 

How TO Compute Critical Exponents. 

The correlation length, is rigorously defined through the correlation function (in quantum 
field theory, the two-points Green's function). For the present discussion it is enough to consider 
£, as the size of the biggest structures in the system in units of the lattice spacing (effective UV- 
cutofF). Thus, if the RG transformation really yields the long scales physics invariant, it must be 
true that 






whichever the prescription may be, as long as the blocking parameter is 2. At the fixed points the 
correlation length may only be zero (complete order or complete disorder) or infinite (nontrivial 
scale invariance, fractal structures). 
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The correlation length is infinite only for the Q = Qc case. Its divergence as we approach that 
critical point fulfills a power law: 

^^(g-gc)-" [5] 

Let us consider a point g quite near g^, so we may consider the R(g) function approximable 
by a straight line with slope R'(gc). Defining the reduced coupling constant g = (g — gc), we may 
say that R(g) ^ R'(gc)g if g is small. 

Thus, when £, £,/2, the coupling constant g — > R'(gc)g- Merging both relations with the 
expected form [5] we obtain: 

lTl2 

^" InR'(gc) 

These computations yield 'v ~ 1.482, while the exact value is just 1. 

Other critical exponents may be computed in analogous ways. We have calculated "v merely 
for illustrating purposes, and the full computations would lead us too much astray. We refer the 
interested readers to [GMSV 95] and the rest of the general works on RG, such as [P JP 82] and [JUL 
81]. 

The conclusion to be drawn, in any case, is that, despite the qualitative general agreement, 
critical exponents are not as precise as it was expected. 



2.2. A Simple Test Model: the Particle in a Box 

As it was stated in the introduction, despite the theoretical appealing of the RSRG techniques 
and their impressive success on impurity problems, their numerical application was slightly disap- 
pointing and was discarded as a practical numerical technique during the 80's. An informal talk 
by Kenneth G. Wilson in 1986 has become legendary [WHI99]. In it he proposed to focus on a 
simple problem which seemed to contain the core difficulties: a quantum spinless particle moving 
freely in a box. 

Once a physical scenario has been settled, a set of robust observables should be specified as an 
appropriate target. In our case we shall focus on the energy spectrum and its associated scaling 
exponents which shall be defined in due time. But in order to gain a deep physical understanding 
(which is the true aim of this work) we shall attempt to depict the full wave-functions. In fact, 
very interesting conclusions shall be drawn from geometrical considerations on these pictures. 

The mathematical formulation of this problem is rather simple. The configuration space where 
the particle dwells is discretized into a graph, where vertices represent cells in real space and edges 
take into account the topology (i.e.: connectivity). The Hilbert space gets finite-dimensional and, 
if no potential is present, the hamiltonian shall be proportional to the laplacian on the graph^. 

The boundary conditions (b.c.) which are imposed on the problem shall prove to be of outmost 
importance. If fixed b.c. are held (U'CO) = A>0-) = 0)) then the hamiltonian is proportional to the 
matrix: 



For further explanations the reader is referred to appendix A. 
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Hf = rl - A 

where r is the connectivity of the bulk vertices, I is the identity matrix and A is the adjacency 
matrix of the graph (i.e.: element (i, j) takes value +1 when vertex i is linked to vertex j and zero 
otherwise). For example, for a one-dimensional finite box split into N cells we have: 



/ 2 



Hf 



V 



-1 
2 / 



This matrix may be exactly diagonalized. The eigenvalues and eigenvectors are given by: 



"'^"^'^^'(^^Tiy' ^ith3e[l...N], nG[0...(N-l)] 



On the other hand, if the boundary conditions are free = M^'CL) = 0), the hamiltonian 

is directly proportional to the combinatorial laplacian on the graph^: 



Hi 



/ 1 -1 
-1 2 



V 



2 -1 
-1 1 / 



This operator has a zero mode: a spatially uniform wave-function is eigenstate of the H; 
matrix with eigenvalue Eq = 0. The rest of the spectrum may be exactly obtained too: 



cc cos (^0 - 1/2)) En = 4 sin^ (^) with j G [1 . . . N], n G [0 . . . (N - 1)] 

There are many different physical problems which lead to the same mathematical formulation. 

It is a good technique to have in mind the alternative physical realizations of a given mathematical 
model when trying to solve it. In our case, a vibrating string or a tightly bound electron in a 
lattice of similar atoms are some of the available possibilities. 

Such Tight Binding Model (TBM) is a consistent analogy. Let us consider a molecule or solid 
(i.e.: the lattice), composed of either atoms or (small) molecules (i.e.: the sites). An electron may 
only occupy a single orbital per site, whose energies are the diagonal elements of the hamiltonian 
or self-energies. But orbitals at nearby sites overlap, thus providing hopping terms (non-diagonal 
elements). Because of them, the ground state of the electron gets dispersed throughout the lattice^. 

Of course, a difFerential condition may not be imposed directly on a discrete structure. See appendix A for 
the accurate meaning of the analogy. 

^ Unless some sites might have a rather low self-energy and act as deep wells, as it shall be explored in section 
4.7. 
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It is remarkable the exposition of [FEY 65], where the TBM is used to introduce the wave-function 
concept. 

BRG FOR THE Particle in a Box. 

Let us consider a ID lattice of N sites: Li = [1 . . . N] and a new lattice which stands at the 
right of this one: L2 = [(N + 1), . . . , (2N)]. The ground states for each of these lattices are known. 
Now we are asked whether we might variationally build up the ground state of the composite lattice 
L = [1 . . . 2N] using the small block ground states as an Ansatz. 

According to the standard BRG technique, we should write an effective hamiltonian for the 
subspace spanned by the two subsystem ground states, extended appropriately to the whole lattice 
by writing zeroes outside its domain. 

The problem reduces to the diagonalization of a 2 x 2 effective hamiltonian matrix. Let \^\>l) 
be the left ground state and \^\>■R) the right one. Thus, Ht shall denote the total hamiltonian (for 
the full system): 

/(tpL|Ht|ll)L) (tpL|Ht|ll;R)\ 
V(^pR|Ht|ll>L) (t^R|Ht|tl>R)y' 

Variational approaches of physical problems are always highly dependent on the quality of the 
Ansatz. In this case, it proves to be surprisingly inappropriate. 

Checking the numbers for a definite size of the system (Nt = 40 sites for the composite system): 



Tl(0) 
■^20 


(ll^LlHt |tpR) 


Tl(0) 
■^20+20 


i:(0) 


0.0223 


-0.00212 


0.0202 


0.00587 



Table 1. Results obtained with the BRG for the free particle in a box. E^j^' is the exact ground state energy 
for a system of N sites, and E^j^^l^j^^ means the approximate energy using a BRG algorithm with blocks of sizes 
Ni and N2. 

The best way to show the cause of the failure of the BRG (400%!) is to plot the wave-functions 
of the exact ground states for the small blocks and the resulting approximation, as it is done in 
figure 5. 

The boundary conditions force the wave-functions to take the value zero at the borders of each 
block, thus making a spurious "kink" appear in the center of the complete system. 

Therefore, the first lesson to be obtained is that "Boundary conditions may be determinative 
for the failure of a RG-scheme^' 

New RSRG techniques were developed which dealt with these problems successfully by the early 
90's. These works pointed out the correlation between the blocks as the reason of the failure of the 
BRG and designed a thoroughly new RG algorithm, called the Density Matrix Renormalization 
Group (DMRG)^. A more insightful analysis of the failure of BRG appeared in 1995 [MDS95], 
which proved its correctness in the design of a modification which, while conservating the old block 
constructions idea, yielded the energies precisely. 

^ It is interesting to remark that S.R. White and R.M. Noack, creators of the DMRG, had formerly invented the 
"combination of boundary conditions", which used as bricks to build the ground state of the fixed b.c. laplacian 
the states of small blocks created with both types of b.c. The algorithm worked, but the method was rejected 
since it was not generalizable to interacting systems [WHI gg]. 



26 2. The Correlated Blocks Renormalization Group. 




Figure 5. Above: the exact ground states of the small boxes are depicted. Below: the BRG prediction for the 
ground state of the full system. 



The original approach used free boundary conditions because of the homogeneity of the ground 
state (which is flat). This fact makes the RG prescription recover exactly the ground state of the 
composite block. The energy of the ground state and the first excited state are obtained with 
the BRG technique with a good precision. At the end of the process, the wave-functions were 
reconstructed and depicted. 

The techniques which appeared in [MDS 95] were further developed and expanded in our work, 
which shall be exposed in the second part of this chapter. 



2.3. Correlated Blocks in ID. 
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Part II. The Correlated Blocks Renormalization Group. 



Among the successful RSRG techniques, the one known as Correlated Blocks Renormalization 
Group (CBRG) is probably the one which resembles the old BRG more closely. The seeds which 
grew to make up the technique were published in [MDS 95] and became a full fledged algorithm in 
[MRS 96], as a part of the present work. 

Along this second part of the chapter, the simple problem exposed at the end of the previous 
section is addressed with a technique which derivates from old BRG, but which takes the correlation 
between blocks explicitly into account. Despite the simplicity of the problem under study, it is 
possible to obtain conclusions which are generalizable to more complicated ones. 

The essential part of the material of this section is published in [MRS 96], but many results 
and discussions appear in this thesis for the first time. 



2.3. Correlated Blocks in ID. 

Let us consider a linear chain of N sites. This chain is split into small blocks of m sites, 
denoted by bp with p G [1 . . . N^] (Nt, = N/m), as depicted in figure 6. 



Figure 6. Splitting a 9-sites chain into three blocks (m = 3). Notice that some links are missing. 



The hamiltonian for each of the fully isolated blocks shall be denoted by Ap. Due to their 
isolation, they shall always have free boundary conditions at both extremes. These matrices shall 
be denoted as self-energy operators. 

In order to re-compose the full chain, the whole set {ApjJ^'' is not enough. The resulting 
hamiltonian matrix would then be (setting m = 3 for definiteness) : 



A, 





/ 1 -1 






-1 2 -1 








0-11 










1 -1 








-1 2 -1 


An J 






-1 1 




V 





/ 

It is therefore necessary to complete the Ap matrices with new kinds of matrices which are 
related to the missing links. Let us remark that the needed completion amounts to (-l-l)'s added 
along the diagonal in the block border vertices and some non-diagonal (— l)'s, corresponding to 
the removed links. 
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Denoting by (pR) = (p p + 1) the (directed) edge linking the block p to p + 1 (or, equiva- 
lently, (pL) = (p p — 1)) we may define: 




if p < Nb 



if p > 1 



which shall be called influence matrices, since the blocks receive through them influence from their 
neighbours. The physical meaning of the ^{pR) matrix is the energetic perturbation which the 
block p receives from the block which stands at its right side. Standard BRG, therefore, would 
have set the intra-block hamiltonian to 

(Hb)p = Ap + B^pL) + B^pR) 

for bulk blocks p G [2 ... (N — 1)]. But the non-diagonal elements corresponding to the missing 
links are still out of our matrix. We introduce the matrices C^p p+i^ for (undirected) edges: 



C 



{p,p+i> 




These matrices shall be called interaction matrices, because they are not to be added to the block 
hamiltonians, but to stay between two blocks. The full hamiltonian may be reconstructed in this 
way: 





/Ai + B(iR) 


C<1,2> 








A2 + B^2L> + B(2R) 


C(2,3) 


Hi = 




^(2,3) 


A3 + B^3L) + B(3R) 













\ 




C(N„-l,Nb> 


'-(Nb-l,Nb) 


ANb +B(NbL)/ 



The structure of this matrix may be represented pictorially, as it is done in figure 7: 
Ai _ A2 _ A3 _ A4 




Figure 7. Pictorial representation of the decomposition of the full free b.c. hamiltonian. 
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The CBRG step proceeds now by forming superblocks by joining the blocks into pairs. A 
superblock hamiltonian is specified by the matrix: 

^Sb ^/'^P+B(PR> C(p,p+i) \ 

The key point is to "renormalize" it into a single block with its corresponding A matrix. This is 
graphically represented in figure 8: 




B(1R) 



Figure 8. The superblock is formed with the elements depicted in the left graph. It shall be renormalized into 
a single new block along with its self-energy matrix. 



As m is the number of states for each block, the superblock has dimension 2m x 2m. That 
matrix is exactly diagonalized, and the m lowest energy eigenstates |i|>i)... |i^m) are retained and 
written down as the rows of a matrix T, which is called the truncation operator: 

m / \ 

T[p,p+i] = = ; 

where |5i) means the m components vector which has a single 1 at position i-th. 

Therefore, T[p p+i] is an operator which takes from the direct sum of the vector spaces of both 
blocks to the m-dimensional space where each dimension "represents" a different eigenstate of the 
superblock hamiltonian. Its matricial representation, thus, has dimension (2m) x m. Its action, in 
physical terms, is to take a vector which dwells on both blocks and return its m weights on each 
of the retained states. 

Its adjoint operator, Tj^ p+i] ^® called the embedding operator, which takes a state expressed 
by its weights on the set of retained states (m numbers) and returns its full form (2m numbers). 
We have: 

m 

t=l 



Now the superblock shall be renormalized to be a simple block in the following RG step. 
Therefore, it shall have an own self-energy operator A^, with p' = which shall be defined as 

the reduction of H^'' to the linear space spanned by the retained states: 



Ap, = THpp^,p_^i]T^ 



30 



2. The Correlated Blocks Renormalization Group. 



So, in components: 



a;, = |6i)(^^i|HS,%+,]|^^^)(6j 



This transformatin is just a "lossy" basis change. Performing this procedme with aU the Nb/2 
superblocks we get a new set of A matrices, which shah be denoted by {A^^}^^'^^. 

The power of the RG stems from its iterativity, i.e.: the possibihty to proceed recursively. In 
our case, we miss a set of renormalized influence (B) and interaction (C) operators in order to 
proceed. 

The renormalized influence and interaction operators are obtained by building a (4m x 4m) 
inter-superblocks eflFective hamiltonian, which contains the links which were missed when building 
the superblock: 



H 



Sb-Sb 
[P...P+3] 





B((p + 1)R> 


> 

C(p + l,p+2) 


*^(p + l,p+2> 


' 



\ 



It is intuitively clear that the B and C matrices should be extracted from this inter-superbloeks 
hamiltonian. But it is necessary to take care in this process, since we have two different truncation 
matrices. Thus, we define the 2m x 4m combined truncation matrices: 



T[p...p+3] = T[p,p+i] 



Ir 



[p+2,p+3] + I[p,p + 1] W l[p+2,p+3] 



I.e.: they act on a field with values on the sites of the four blocks and returns 2m numbers. The 
first m are the weights on the left part, and the last m are the weights on the right part. 
The inter-superblocks matrix is truncated using these combined operators: 



L|Sb-Sb yt _ / (P'R) "^(P'.P' + I) 

'[p...p+3] ri[p...p+3] '[P...P+3] - 1 p't ^f 

So we may read the new sets of matrices: 

l°{pL>/2 l°(pR)/l 



Nb/2-l 



j-w -,Nb/2-l 



From this point on, the CBRG iteration is easily established. The main idea is to recur the 
procedure until the whole system is of the size of a superblock, rendering therefore the approximate 
diagonalization of the full hamiltonian feasible. 

The details for a concrete numerical implementation are left for the next paragraph. 



2.4. The CBRG Algorithm in Practice. 
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2.4. The CBRG Algorithm in Practice. 

Throughout this section, the technical details of the main applications of the CBRG are dis- 
cussed. 



Implementation of the Homogeneous System. 

It is usual to discard the computational aspects when discussing a mathematical-physics al- 
gorithm. It must be remarked that without the concrete application within a computer language, 
the above discussion is almost void. For a general discussion on the numerical implementations see 
appendix B. 

Since the system is homogeneous, all the A| matrices are equal, and the same holds for the 
^(iL)) ^(iR) ^-i^d C^p q). Therefore, the data structure is rather simple: the system just consists 
of a "pack" {A, B]_, Br, C}. Let us remind that m is the number of sites of the original blocks 
and Nb = N/m is the number of such blocks. Therefore, matrices A, Bi_, Br and C shall have 
dimension m x m. 

A procedure called Composition is of outmost importance. It builds a 2n x 2n matrix out of 
four n X TL matrices. It works like this: 

. . / ^ Ai A2 

Composition (A1,A2, A3, A4) = 

\ A3 A4 

The superblock hamiltonian 2m x 2m is written as (Ct means the transpose of C) : 

(A -f- B C \ 

This superblock hamiltonian is diagonalized. The m lowest eigenstates are retained and the 
truncation matrix T (m x 2m) is written down with these vectors as rows. 

The easiest renormalization procedure is that for the A matrix. If Tt denotes the transpose of 
T (i.e.: the embedding operator), then^: 

a' ^ T * Hsb * Tt [6] 
The rest of the matrices renormalize in this way (Z denotes a null m x m matrix): 

HOO = Composition (z, Z, Z, BR) HOI = Composition (z, Z, C, z) 
HIO = Composition (z, Ct, Z, z) Hll = Composition (bl, Z, Z, z) 

BR' ^ T * HOO * Tt 

C' <- T * HOI * Tt [6'] 
bl' ^ T * Hll * Tt 

(Of course, HIO is not needed since it would renormalize to Ct, but it has been written down 
so as the pattern is more easily recognized). 



^ We shall employ the notation x ^ A to denote the assignment of the value A to the variable x, as it is common 
practice in computational sciences. 
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The number of superblocks to be diagonalized is just the number of RG-steps which must be 
taken before the superblock contains the whole system, i.e.: log2(Nb/2). 

Numerical Results for the Free Homogeneous System. 

In the homogeneous chain with free boundary conditions the results are satisfactory. Table 2 
shows the numbers for a chain with free b.c. at both ends with 768 sites split into blocks of 6 sites 
each (Therefore, N^, = 128 and 6 RG-steps are required). 



Energy 


Exact 




CBRG 




Eo 







1.1340 X 10- 


14 


El 


1.6733 X 10- 


-5 


1.9752 X 10- 


-5 


E2 


6.6932 X 10" 


-5 


7.6552 X 10" 


-5 


Es 


1.5060 X 10" 


-4 


1.8041 X 10" 


-5 


E4 


2.6772 X 10- 


-4 


2.9681 X 10" 


-4 


E5 


4.1831 X 10" 


-4 


5.1078 X 10" 


-4 



Table 2. Numerical results for a free chain of 6 x 2^^ = 768 sites. 



The results for the excited states, with errors in the range 15 — 25%, are quite correct quanti- 
tatively. We shall obtain more robust tests, such as scaling exponents of the energies. 

The ground state energy is always zero, which is reproduced within machine precision. In the 
large N regime, the excited states fulfill a simple relation: 



,2 



00 



A least squares fit of the data to a power law of the type 



yields a = 2.09 ± 0.03 and (3 = 1.95 ± 0.03 using u G [0 . . . 6] and N G [96 . . . 3072]. The biggest 
source of error lies in the prefactor, which fits to K = 7 ± 1 (the exact value is ^ 9.86). 



Inhomogeneous System. 

The data structure for the system when it is inhomogeneous is slightly more intrincate. Every 
Ai, B^ii_^, etc. matrix is different. At each RG-step, memory should be reserved for a different 
number of matrices: matrices of size m x m are stored at the beginning, but this number is 
divided by 2 after each step. 

The only step which is slightly dangerous and should be kept in mind is that the T matrices 
employed in equations [6] and [6'] differ. 

For definiteness, let A(i) be the matrix corresponding to the i-th block at a given RG step, 
while BL(i), BR(i) and C(i,i+1) denote the corresponding influence and interaction matrices. 
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Let be the number of blocks at this RG step. Therefore, Nb/2 superblock hamiltonians should 
be written down. For any i € [1 . . . Nb/2], 

Hsb(i) = Composition (A(2i - l) + BR(2i - l),C(2i- l,2i),Ct(2i- 1, 2i), A(2i) + BL(2i)) = 

_ /'^2i-l + B((2i-1)R) C((2i_i),2i> \ 
V C|(2i_i),2i) A2i + B^(2i)L>y 

Diagonalizing this superblock hamiltonian and discarding the highest energy eigenvectors we 
obtain a truncation matrix T(i). The renormalization of the A matrices is straightforward: 

A'(i) <— T(i) * Hsb(i) * Tt(i) 
The inter-superblocks hamiltonian is written this way: 

BR'(i) ^ T(i) * Composition (Z, Z, Z, BR(2i)) * Tt(i) 
BL'(i) ^ T(i + l) * Composition (BL(2i + l), Z, Z, z) * Tt(i + l) 
C'(i, i + l) ^ T(i) * Composition (z, Z, C(2i, 21 + l), z) * Tt(i + l) 

The number of RG steps given before the superblock hamiltonian covers the whole range is 
stilF log2(Nb/2). The number of superblocks to be diagonalized is now given by: 

logj Nb/2 

1 + 2 + 4 + - •• + Nb/2 = 2^ = 2Nb-l 

t=0 

According to [PTVF97], typical diagonalization of a tridiagonal matrix takes O(rL^) steps 
(without eigenvectors). Our technique, which is totally general for tridiagonal matrices, gives a 
good numerical estimate in only 0(n) steps (0(log2(n)) if the homogeneous algorithm may be 
applied) . 

Fixed boundary conditions. 

The case of a homogeneous chain with fixed (or mixed fixed-free) boundary conditions may be 
included under the title of inhomogeneous system. As it is discussed in appendix A, the natural 
system is the free chain (isolated system) , while the fixed chain implies a trivial yet existent exterior 
space, which makes sites near to the border "really" different. 

This issue is made more clear in our approach. In order to solve the fixed boundary conditions 
system it is necessary to include B^^l^ and B^jsi^^r^ matrices, which account for the influence of the 
(trivial) sites outside the chain (the so-called "tack" sites). 

Of course, it is not required to compute all the superblock matrices. The system is almost 
homogeneous and therefore it is enough to treat differently the border blocks. 



^ If a RG-step index is added to the matrices (e.g. A^'^^(i), C^'^^(i,i + 1)...) and all the levels are stored, the 
resulting data structure may remind strongly of wavelets [lem 8g]. Of course, this fact is not casual: see 1.4. III. 
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Potential Error range 



V(x) 


= with free b.c. 


0- 


22 % 


V(x) 


= with fixed b.c. 


17- 


-28 % 


V(x) 


= with fixed-free b.c. 


17- 


-22 % 


V(x) 


= sin(8x) + cos(84>x) 


0.02 - 


- 0.2 % 


V(x) 


= 1/2 (x- 0.5)2 


3- 


-5% 


V(x) 


G 7^[o,l] 


10-4 


0.4 % 



Table 3. Some results for inhomogeneous ID systems, with 768 sites and 6 states. In all cases, the variable 
X = {i — Vj/Nt, i.e.: the site index over the number of sites, (p = 0.5(1 + \/5) is the golden section and V € 7?,[0, 1] 
means that V{x) is a random variable equally distributed on [0, 1], drawn independently for each x. 

Numerical results for Inhomogeneous Systems. 

Any potential may be used to test the CBRG in a non-homogeneous ID system. Many of 
them have been checked, and a few results are shown in table 3. 

A more robust check shall be performed on the fixed boundary conditions and the harmonic 
oscillator potential V(x) = 1/2x2. ^.j^g gj,g^ case, a (u + 1)2/N2 law may be obtained. The 
second is characterized by the equally spaced spectrum. 

• Fixed boundary conditions. 

Fitting to a functional form E]^ = K(n + 1)'^/NP we obtain K = 7 ± 1, a = 2.05 ± 0.02 and 
P = 1.94 zb 0.02. The N and m values are taken from the same range as in the free case, and the 
results may be seen to be similar in both cases. The exact values in the asymptotic regime are, of 
course, a = (3 = 2 and K = 7t2. 

• Harmonic oscillator. 

We fit the energies to the law = K(n + l/2)'^/Nl^. The asymptotic exact values are, in 
this case, a = (3 = 1. The numerical calculations for the same range (6 retained states and size of 
the system from 192 to 3072) yields a = 1.011 ± 0.008 and |3 = 0.988 ± 0.007. 

Many questions related to the relative success of CBRG in so different systems need an answer. 
We shall postpone the due explanations until the wave-functions have been depicted. 

Fixed points of the CBRG. 

Let us consider the ID homogeneous case with free b.c. The parameter space for the RG 
procedure is formed in this case by the matrix elements of A, Bl, Br and C. Which are the fixed 
points for these parameters? 

Of course, there is a trivial fixed point corresponding to old BRG: Bl = Br = C = 0. 
The evolution of matrix A towards this fixed point is interesting. Each RG-step the m/2 lowest 
eigenvalues of matrix A are duplicated and the rest is discarded. Thus, matrix A reaches the null 
value in a finite number of steps even though its eigenvalues never decrease in magnitude. 

The non-trivial fixed point is more interesting. Let us take, e.g., m = 3. Then all matrices A, 
Bl, Br and C tend to zero asymptotically according to power law: 
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A*=N-^a*, Bt = N-^bt, = N'^b^, C* = N" 
And the matrices in the fixed point are 



a 












(: 


ei 


















— s 






— s 


t 








t 








With ei « 2.3, 62 ~ 9, s 1.3 and t ^ 1.8. The fit for z and z' yields numbers which are 
compatible with 1 and 2 respectively. Thus, starting with the matrices a*,... and applying 2^ or 
2^ after each RG step, they stay the same. 



2.3. Wave-function Reconstruction. 

One of the most interesting features of the CBRG technique is the possibility of wave-function 
reconstruction. It is not only an interesting property by itself, but it shall also lead us to a deeper 
insight into its mechanisms. 

The key to the procedure is the storage of the truncation operators at each step for all su- 
perblock hamiltonians. They contain the full information for the reconstruction of the wave- 
function, even if it is stored in an intrincate recursive fashion. 

The first RG step is the only one to take place in real space. All other steps take place on a 
space where only the difference between left and right is real, but where the states for each part 
do not represent real sites. They stand for the weights of the smaller scale wave-functions on its 
side. 



An illustrating example. 

We shall exemplify the procedure with the most simple non-trivial case available. Let us 
diagonalize the block hamiltonian with 4 sites. The two lowest energy eigenvectors are interpreted 
as two real space functions which are shown in figure 9A. 



Figure 9 A. The two lowest energy states for a 4 sites free-free chain. 
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Now we form the truncation operator using the two shown eigenstates as rows. The fohowing 
RG step proceeds by joining two systems similar to the former one. The new system has a left 
block and a right block, but its "sites" only represent the weights of the old states on each of the 
sides. 

Let us form, thus, the effective superblock hamiltonian for the four mentioned states (two for 
the left side and two for the right side). Such a superblock, even though it represents 8 real sites, 
is only a 4 x 4 matrix. Its 2 lowest energy eigenstates are shown in figure 9B. 



New ground state 



New first excited 




Ground Exc. Ground Exc. 




Ground Exc. Exc. 

Ground 



Figure 9B. The two lowest energy eigenstates of the superblock. The left graph represents the formation of the 
new ground state, and the right one that for the new excited state. Bars represent the weights of the old block 
states. 

The depicted states show that only the left and right ground states contribute to the ground 
state of the bigger system. Both states are flat, and so is the global state. 

The first excited state is more interesting. The left side receives a positive contribution from 
the small ground state and another from the small excited state. This implies that it shall have the 
same shape as the small excited state, but ^^raised" due to the constant term. On the right side, 
the contribution from the ground state is negative, which implies a uniform ^^descent" to the other 
contribution, which is from the old right excited state. Figure 9C shows the full wave-functions at 
the end of the last step. 



Figure 9C. Full wave-functions at the end of the RG step. 



2.3. Wave-function Reconstruction. 
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The Reconstruction Algorithm. 

If all embedding matrices are stored, then the complete wave-functions may be recon- 
structed. The T,*^ matrices take a vector of m components (which represent the weight of each of 
the m lowest energy states of the superblock), and return a 2m components vector. The first m 
components of this vector represent the contribution from each of the left states, meanwhile the 
last m represent the contribution from the right ones. 

The reconstruction algorithm works downwards. In a certain sense, it reminds of the appli- 
cation of more and more powerful magnifying glasses. Starting from the eigenvectors of the last 
superblock hamiltonian (which represent the whole system), the previous step is applied to 
them. Since that operator doubles the number of components, we may introduce two operators 
called L and R, which represent respectively the restriction to the left and right halves of the 
system. Matricially, 




Figure 10. This tree shows the reconstruction algorithm. An initial superblock eigenstate taken from the last 
RG-step is subsequently split into left and right parts until the basic level (real space) is reached. 

Figure 10 shows how the reconstruction tree develops. Its root is the highest level and it 
grows until its branches reach ground: real space. If each R is written as a and is L as a 1, the 
prescription to obtain a given block may be easily given for the homogeneous case. 

• Write the block index in base 2. The number must have binary digits (i.e.: bits), the 
highest of which may be zero. 

• Use the dictionary 1 L, R to obtain a "word" in the {L, R} alphabet. 

• Multiply all matrices, taking each L or R from the appropriate level. The columns of the 
resulting matrix are the real space components of the global states on the sites corresponding 
to the given block. 

In the inhomogeneous case the only difference is that at each level there are more than a pair 
of {L, R} matrices. A general formula may be provided. Let ig € [0 ... N — 1] be a site in real space. 
Thus, ii = LVttlJ is the index of the block it belongs to, and jo = i-o ~ Tni-i is its position index 
within that block: jo £ [0 ... m — 1] and li G [0 . . . N^, — 1]. Following this line, we may define a 
whole set of indices: 



38 



2. The Correlated Blocks Renormalization Group. 



. ^ / L^k-i/2j ifk>2 
- \ [io/mj if k = 1 

• = f^k- 2i-k-i if k > 1 
-^"^ ^ \ io - mil if k = 

Thus, i]^ is the index of the k-th level block in which the site is included, and jic is the sub-block 
within that block. At levels > 1, ji< may only take the values and 1, standing for right and left 
respectively. 

The last convention we shall adopt shall be to define B(l,i, 0) to be the right part of the 
embedding operator at the l-th RG step for the i-th block. Consequently, B(l,i, 1) accounts for 
the left part. 

With this notation, if there are levels (m x 2'^'- = N), the io-th component of the l-th 
global wave-function is given by: 

< = [B(0,ii,ji)].^l,jB(l,i2,j2)k,,^,---[B(N;-l,iN,-l,jN,-l)]„,^_,,; [7] 

where the summation convention is implicit on all the k^ indices. 
Graphical results. 

Once the reconstruction algorithm has been described, we shall reproduce some of the wave- 
functions which have been obtained with this technique. Figure llA shows the most simple case, 
the homogeneous algorithm for a free chain of 64 sites. 




Figure 11 A. In continuous line, the CBRG wave-functions for the free system with 64 sites. The dots mark the 
exact solutions. 

The ground state of the system is exactly reproduced, as it might have been foreseen. The 
excited states may be seen to "degrade gracefully" as the energy grows. The number of nodes is 
conserved, which is the key for the accuracy of the results. It is always true that states with even 
symmetry stand a better reproduction (and, therefore, smaller errors in energy), because the zero 
slope is respected at the origin. 



2.3. Wave-function Reconstruction. 
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Figure IIB. The fixed b.c. wave-functions with the same plotting conditions as in figure llA. 




Figure f 0. Lowest energy states for the harmonic osciUator according to CBRG. 

Figure IIB shows the full wavefunctions for the system with fixed boundary conditions (without 
potential). Figure IIC shows the lowest energy states for a system with free b.c. but in presence 
of a harmonic potential V(x) = (x — 1/2)^ with x G [0, 1]. 
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2.6. Two-dimensional CBRG. 



The CBRG algorithm is by no means hmited to work in one dimension. Influence and interac- 
tion matrices are easily generalized to work in higher dimensional systems. In this work we shall 
only expose the 2D case. 

Let us consider a 2 x 2 block with free boundary conditions. The laplacian for such a graph is: 




/ 2 
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2 





2 



V 
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-1 

-1 
2 / 



Influence and interaction are related to neighbouring blocks, which are four in our case. From 
the graph theory viewpoint, influence B matrices are associated to directed links, and interaction 
C matrices to undirected links. 

Each block must now have four C matrices (CdU) CuD) Clr and Crl) and four B matrices 
(Bu, Bd, Br and Bl), whose role is shown in the the following graph, along with their numerical 
values: 
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It is easy to convince oneself that Clr = C^j_ and Cdu = ^uu- Under fixed boundary 
conditions, the laplacian would have been: 



Lfixed = A + Bl + Br + Bu + Bd 
The superblock is formed by four blocks, and its hamiltonian is a 16 x 16 matrix in our case 



/ 4 -1-1 \ 

-14 0-1 

-10 4-1 

V -1-1 4 / 
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(in general, 4m x 4m). Its graphical representation coincides with the graph for the small 2x2 
system so its matrix expression is: 



H 



Sb 



/A + Br + Bd 
Crl 
Cdu 

V 



Clr 
A + Bl + Bd 


Cdu 



Cud 


A + Br + Bu 
Crl 



\ 

Cud 
Clr 
A + Bl + Bu/ 



The superblock hamiltonian is diagonalized and the truncation and embedding operators, T 
and T'^, are obtained using the same operation as before. With these 4 x 16 and 16 x 4 matrices 
we start the renormalization of all operators: 

A' = THsbTt 

In order to renormalize the rest of them we need to write down the inter-superblocks hamil- 
tonian. This should contain 2x2 superblocks, but it is easier to write it down in parts: a vertical 
and a horizontal inter-superblock. The first one links two superblocks on a column: 
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The horizontal inter-superblocks hamiltonian is calculated in the same way: 



H, 



hor 

Sb-Sb 



II hor 
•^Sb-Sb 
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From these 32 x 32 matrices we may read the new influence and interaction operators. If we 
denote, e.g., H'^°^(DL) the down-left quarter of the H5|^_[5|^ matrix, we have: 



B'o = TH^^^*(UL)Tt C[ 



UD 



DU 

b; 



c 



RL 



-^vertpgjt 

= TH^°^(UL)Tt 
-- TH^°^(DL)Tt 



b; 



u 



LR 



B'l 



= TH^°^(UR)Tt 
TH^°^(DR)Tt 



This way, the RG cycle has been closed. 
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Wave-Function Reconstrunction. 

The reconstruction of the wave-function is also possible in the 2D case, but it is certainly more 
involved. The global idea is the same, but the details are more cumbersome. 

An auxiliary idea which may help to understand the procedure is that of a quad-tree. It is a 
usual concept in fractal image compression [FIS 95] and consists in a data structure which is able 
to store a 2D image in such a way that: 

• Information is local. 

• It is easy to establish a UV-cutoff in the image (so as to compress it). 

• This UV-cutoff may depend on position. 

The first condition is important so as to obtain high compression rates, but usual Fourier- 
transform algorithms (such as JPEG), storing only global information, suflfer from Gibb's phe- 
nomenon (wave-like distorsion near edges) and other undesirable effects: aliasing^. 

Figure 12 shows a typical quad-tree structure spanning a square region. The outmost node is 
split into four links, each of which is also split into four minor ones and so on. 




Figure 12. A graphical representation of a quad-tree structure. 

If nodes in a block are numbered as in the first picture of this section, the quad-tree structure 
allows us to substitute the binary notation used in the ID algorithm with a base 4 notation. Each 
block is determined by a word in the alphabet {0, 1, 2, 3}, where each letter is to be read according 
to the following dictionary: 

<r-^ Upper-Left 1 Upper-Right. 
2 Lower-Left 3 ^ Lower-Right. 

The embedding matrices may be split into four matrices, corresponding to upper-left, upper- 
right, lower-left and lower-right regions of the superblock. Let them receive names from Bic(O) to 
Bk(3) (with k denoting the RG-step). An equation analogous to [7] may be written down for this 
case. 

Let X and y be the coordinates of a site in real space, Xj-y G [0 . . . L — 1]. Let us also consider 
the smallest blocks to be of size 1^ x l^,, and the number of RG-steps to be N;. We define the 
following sets of integer indices: 



The term aliasing refers to the set of distortions which are introduced in a digital image when some of its 
spatial frequencies with a high weight are (almost) commensurable with the sampling frequency. 
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, _ r [x/lbj if i = 
''^"U<-i/2j ifi>0 

I.e.: the x-coordinate of the i-th level block to which the site belongs. The same definition goes 
for the iji. These and indices are only needed so as to define the following: 

_ J X — Xolb if i = 
''^ " \ xt - 2x|_i if i > 

This is to say: the position inside each i-th level block of the lower level block to which the site 
belongs. The yt are defined in an analogous way. Now the set of indices jt for i, > is defined: 



k = 2yi + Xi 

These jt indices specify the quad-tree path to find the site from the root node. In effect, the 
analogue of [7] is: 



Numerical Results of the 2D Algorithm. 

The procedure described in the previous section is the most simple 2D-CBRG algorithm pos- 
sible. It can be easily generalized to include a potential or non-trivial boundary conditions of any 
kind. 

For the 16 — > 4 RG scheme proposed above, energy stays within 10% of the exact value for all 
states in 6 RG-steps (128 x 128 sites). The scaling relation is almost perfect: the energy fulfills 
E]^ oc N"^ with a ^ 1.99999981. 

Figure 13A shows the four lowest energy states for the 128 x 128 lattice. The third excited 
state is also represented in figure 13B within an axonometric projection with contour lines depicted 
on the basis plane. 




Figure 12 A. Density plot for the lowest energy states for the free b.c. laplacian on a square lattice of 128 x 128 
sites. 
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Figure 13B. Axonometric plot of the third excited state for the free b.c. laplacian on the 128 x 128 lattice, 
obtained with CBRG. 



This section might also have been entitled "Why do free boundary conditions work for ho- 
mogeneous CBRG while the fixed ones do not?" The answer shall provide us with an interesting 
insight on the full CBRG procedure. 

The correct point of view in order to find the answer is the "bricks" image. Let us focus 
again on the homogeneous ID case. Each RG step works in the following manner: ^^huild the best 
approximation to the lowest energy eigenfunctions for the whole system using as bricks the lowest 
energy eigenfunctions for the left and right halves" . 

Self-replicability. 

We might have started the research the other way round. Given a function £ C[0, a], let us 
define the operators 1-[o,q] and R[o,a] as: 



I.e.: they provide us with reduced copies which are similar to the original function for each of the 
parts (left and right). 

Starting from a function cj) we obtain a pair of them, LcJ) and Rcj). We shall try to reproduce 
the original cf) within our subspace spanned by Lcj) and Rcf). Is this possible? 

Let us consider all functions to be normalized and let us define the replica transformation: 



which is sure to be the best approximation within this subspace. Its accuracy shall be given by 
the parameter 



2.7. Self-Replicability and the CBRG. 





7^4) = 4) 



approx — 



(ct^lLcf)) L4) + (ct^lRcf)) Ret) 
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S = (ctJapproxIcfj) 

where the symbol S stands for self-replicability. The value 1 means "perfect" . 

In practice, of course, we are working with discrete functions stored in a computer. If the new 
functions must have the same number of sites (i.e.: they belong to the same discrete functional 
space), then two values of the old function must enter a single site for each of the parts (left and 
right). The most symmetric solution is to take the average of both values. Thus, 

(L4)), = I 5(*2i-i + ct)2i) ifi<N/2 
otherwise 

along with an equivalent formula for the right side. This is the computational expression we shall 
assume. All the examples shall use this implementation on a ID lattice with 256 sites. 

Numerical Experiments in ID. 

Let us apply the process on the ground state of the laplacian with fixed b.c. on a ID lattice. 
The best approximation is given in figure 14A. 




Figure 14A. The ground state of a ID laplacian with fixed b.c. is not self-replicable. 

Here S takes the value 0.84883, which does not seem to be too low when one observes the 
obvious differences. But the replica procedure may be iterated. Figure 14B shows us the second, 
third and fifteenth iterations. 




Figure 14B. The procedure is iterated. The last box represents the fifteenth iteration. 

After some iterations the finite resolution of the computer yields a quasi-constant function as 
approximation. This function is exactly self-replicable. 

The key idea is that, when a function is self-replicable, it is exactly attainable by a homogeneous 
CBRG procedure (which only requires ^ log2(N) operations!). 
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The procedure is easily extended to sets of functions. Let us denote any such set, with m 
functions, as {c^ili!^!- These functions must be approximated within the subspace spanned by 
the 2m functions {L4)i, R4)i}™ All such functions may contribute to the reproduction of their 
"sisters" . 

For example, the low energy spectrum of the ID laplacian with free b.c. is not exactly self- 
replicable, but it is to a good approximation, as it is shown in figure 15A. 




The S parameters are 1, 0.999613, 1 and 0.995715. This means that the ground state (flat) 
and the third states are exactly reproduced. The weights show that: 

• The first state is absolutely self-replicable by itself. 

• The second consist of two copies of itself, the first one raised (using the first state) and the 
second one lowered (also using the first state). The finite slope at the origin is not correctly 
represented. 

• The third state only requires the second one. 

• The fourth state is even more interesting. Both the left and the right parts are a combination 
of the second and third states. The finite slope at the origin is again incorrectly represented. 

The procedure may be easily iterated without excessive distortion. The results are shown in 
figure 15B. These functions have a rough look, but factors S are not too different: 1, 0.99958, 
0.999589 and 0.99492. 

It is perhaps more illuminating to perform the same experiment on the four lowest energy 
states of the fixed b.c. laplacian (see figures 16A and 16B). 

The values of the S parameters at the first iteration are not excessively bad^: 0.953725, 
0.999997, 0.945245 and 0.99998. But the fixed point yields very different numbers: 0.843138, 
0.872436, 0.851632 and 0.80208. 

The aspect of the fixed points of the replica transformation is often quite rough. Appendix E 
describes the meaning and calculation method for a "quasi-fractal" dimension which corresponds 
to the enery scaling exponent under RG transformations. The value for smooth functions is approx- 

9 2 

Sometimes functions which are quite different may be separated by a short distance in £ .ft might be more 
appropriate to use a discrete Sobolev space (see section 5.4), in order to tafce derivatives into account. This way, 
an overimposed roughness might be detected. 
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imately 2. Values around zero might be physically interesting. The fixed points for the eigenstates 
of the free and fixed b.c. laplacian yield the results of table 4. 



48 



2. The Correlated Blocks Renormalization Group. 



Ground state 1st exc. 2nd Exc. 3rd Exc. 



Fixed b.c. 


0.034 ±0.027 0.21 ±0.08 


0.33 ±0.08 


0.42 ±0.11 


Free b.c. 


— 1.55 ±0.06 


1.65 ± 0.05 


1.4 ±0.05 



Table 4. Scaling dimensions of the energy under RG for the fixed point reached with free and fixed b.c. wave- 
functions. The free ground state is missing since its calculation is meaningless: its energy is exactly zero. 

Generalization to 2D. 

The process may be easily generalized to 2D, if instead of splitting the interval into two regions 
we part it into four. The case of the eigenfunctions of the free b.c. laplacian yields a fixed point 
which is much smoother than in the ID case, as it is shown in figure 17. 




In the fixed b.c. case, we obtain a rather different fixed point, as it is shown in figure 18. 




Figure 18. Same concept as in figure 17, but for initial wave-functions with fixed b.c. The fixed point is 
thoroughly different. 



2.8. Bibliography. 
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Figure 19. Wave-functions 16 x 16 of higher energy than the previous ones, which yield a repeated pattern 
based on the Pascal or Sierpihski triangle. It must be remarked that this curious structure is not stable, and that 
it is due to a slight asymmetry of numerical origin in the initial states. 

Among the fixed points we have found a great richness of structures. Figure 19 shows a pattern 
which may result famihar. 

Self-replicability and the CBRG. 

Once the question around the free b.c. has been settled (they work because they are approx- 
imately self-replicable), there appears a second question: why are free b.c. states approximately 
self-replicable? The answer is: because they resemble a fixed point of smooth functions, towards 
which they do not converge. 

There is a great number of fixed points for this replica transformation-^", but most of them are 
non-smooth. The only such set to the knowledge of the author is really simple: the polynomials. 

It may be proved that the set {l,x,x^ . . .x^} is exactly self-replicable for any n and inter- 
val. The reason is that translations and scaling transformations of the set return as garbage a 
contribution from lower order polynomials. 

Therefore, in a certain sense, we might say that the success of free b.c. with CBRG is due to 
the fact that their eigenstates look like the polynomials. Fixed b.c. states also resemble them, but 
they lack the two first ones: 1 and x. This guarantees their failure. 
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Part I. DMRG for Quantum Mechanics. 



The "particle in a box" problem, in which K.G. Wilson had found a synthesis of the difficulties 
of the BRG approach-^, was solved for the first time using RSRG techniques by Steve R. White in 
1992 [WHI 92]. The technique became known as Density Matrix Renormalization Group (DMRG). 
With its appearance on stage, RSRG received in a few years full acknowledgement as a high 
precision numerical method. Since Wilson's solution of the Kondo problem no RSRG method had 
ever achieved comparable results. In a sense, DMRG holds a certain "magical" aura: its success is 
still surprising both to the sporadic practicioner and the specialist. 

Our work on the DMRG has focused on the less intrincate of its applications so as to reach a 
deeper understanding of its inner workings: quantum mechanics of a single particle. 

The first two sections of this chapter contain a detailed introduction to DMRG applied to 
quantum mechanics. The basic aspects of the DMRG are taken from the original works of Steven 
R. White and Reinhard M. Noack ([WHI 92], [NWH92], [NWH 93] and [WHI 93]) and from [GMSV 



See section 1.4.V and the first part of the previous chapter. 
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95], which contains its first appearance in a "textbook". The apphcation to quantum mechanical 
problems is based on the work of M.A. Martfn-Delgado, R.M. Noack and G. Sierra [MSN 99]. We 
would like to remark that all the numerical results exposed in this work have been obtained with 
our own programs. 

Afterwards, an extension of the method is developed which was carried out by our group [MRS 
00]. Initially, the DMRG was only applied to unidimensional problems, but it is also possible to 
employ it to analyze systems whose configurational space has a tree-like topology. In the second 
part of this chapter a suitable extension is applied to the study of the dynamics of excitons in a 
kind of polymeric molecules which have been recently synthesized, known as dendrimers. 

The literature on the DMRG is rather wide and quickly growing (see the article by Karen 
Hallberg [HAL 99] for a global point of view). In this thesis we have only cited the works which 
were most closely related to ours, and no importance sampling may be assumed. 



3.1. The Density Matrix Renormalization Group. 

The global idea of DMRG may be stated in a single sentence: "when splitting a system into 
blocks, these should not remain isolated" . 

S.R. White realized that the lowest energy states of the isolated blocks need not be the best 
building bricks to construct the big block. It is necessary for the block to be related to its envi- 
ronment so as it may choose the best states. 

The solution starts with the superhlock idea: a block which is made of some smaller blocks. Its 
dynamics is studied and the lowest energy state is found. From that state, the block states which 
"fit" better may be obtained, i.e.: the best bricks. This objective is accomplished by introducing 
the density matrix. 

The development of this section is freely based on the text of R.P. Feynman [FEY 72] for the 
analysis of density matrices, and the work of S.R. White [WHI 98] for the general DMRG algorithm. 

Why Density Matrices? 

The DMRG algorithm receives its name from one of the tools incorporated into it. It is 
arguable whether the density matrix is the key ingredient of the method, but history has made it 
be considered to be so. 

The reason for which the density matrix appears is the necessity of "fitting" : it is necessary to 
find out which are the block states which reproduce most accurately a chosen global state, which 
shall be termed the target state. The manner in which this state for the superblock is found is a 
problem which shall be addressed later. 

The states which fit better shall be considered to be the most probable ones, in which the block 
may be found, even though they have nothing to do with the low energy states for an isolated block 
(see an example in figure 2). 

Density matrices are, according to some authors, the most fundamental way to describe quan- 
tum mechanical systems. In our description we shall adopt the more "practical" view of R.P. Feyn- 
man [FEY 72], according to which these matrices are necessary because systems must always be 
separated from their environment in order to be analyzed. This comment was the one which took 
S. White to use them as the foundation of his RG formulation. 



3.1. The Density Matrix Renormalization Group. 
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Even though we shall follow [FEY 72], the topic is discussed in similar terms in [LL 71] and in 
[CDL 73]. At the end of the section we shall briefly comment on the other possible view. 

Figure 1 shows the structure of the entire system (the superblock) which is split into block 
and environment. 



Superblock 



Block 



Environment 



Figure 1. Superblock split into block and environment. 



Let the ground state of the full superblock be given by the (pure state!) ket \^\>o)- In a typical 
many-body situation, this ket contains entangled information for block and environment. The 
state is expressed as a sum of tensor products: 

|4'o) = X^M/t,j|ai)®||3j) 

where |ai) make up a full basis for the block and ||3j) one for the environment. Let us consider an 
operator A^ which only acts on the block variables. Then, 

so, since (|3i||3j) = 6ij, we have 

Now we may define the Wocic density matrix through its matrix elements in a complete basis 

(cxt.|p'^|cxt)= J]vl;?.vl;,,,^ 
j 

i.e.: the environmental components are "traced out". With this definition our expectation values 
may be written as: 

The block density matrix is positive-defined, self-adjoint and has unit trace (see [FEY 72]). 
The interpretation, nevertheless, varies with the author to be considered. The original formulation 
was given by L.D. Landau and J. von Neumann [OMN94] in order to introduce thermal averages 
in quantum mechanics. The target was to deal with any ^^second level of uncertainty" beyond 
the one coming from quantum mechanics. Density matrices are a fundamental tool to study 
thermodynamically a quantiim system or to take into account any lack of information concerning 
the exact way in which the state was prepared. In other terms: 

" Whether the environment is not relevant for the computation or if it is out of reach, the block 
(system) must not be considered to be in a pure state, but in a mixed state. I.e.: an statistical 
ensemble of pure states with given probabilities." 
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The states we are referring to do not need to be orthogonal, but each one of them must of 
course correspond to a different state of the environment. We remark Glcason's theorem, which 
asserts that, under very general conditions, every mixed state may be described by a self-adjoint 
positive-defined and unit-trace density matrix [GP 8g]. 

It must be noticed, even though briefly, that the theoretical physicists which have considered 
seriously the measurement problem in quantum mechanics tend to give a different interpretation: 
"A density matrix is the most complete possible description of an ensemble of identically prepared 
physical systems^\ According to this second picture, we should have started our discussion with 
a density matrix for the complete system. For a defence of this position the reader is referred to, 
e.g., [GP89] and [BAL90]. 

DMRG: The Basic Idea. 

The objective of the DMRG is to choose for each block the states which are the best possible 
building blocks. Thus, an eigenstate of the full system is sought, which is called the target state. 
Afterwards we build the set of states of the block which are most appropriate to reproduce that 
state. 

This "fit" is performed by building a density matrix for the block from the target state, by 
"tracing out" the variables which belong to the environment. This density matrix is diagonalized. 
Its eigenvalues are the probabilities for each of the different states of the block. Therefore, it is 
necessary to choose the eigenstates of highest eigenvalues of the density matrix. 

Figure 2 shows graphically the difference between choosing an eigenstate of the density matrix 
(with highest eigenvalue) and an eigenstate of the block hamiltonian (with lowest eigenvalue) . 

Block 




B 



< Superblock > 



Figure 2. A QM example showing how the p eigenstate fits perfectly the target state at the block, while the 
eigenstate, which has to boundary conditions, does not. 

Even though the details of the procedure are strongly dependent on the particular implemen- 
tation, the main lines of work are essentially the ones that follow. The basic DMRG algorithm 
consists of a warmup, which establishes a consistent set of variables, and the sweeping, which 
iterates the process on the full system until convergence is reached. 

The following discussion is rather formal, and all the details on concrete implementations are 
left open. 
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DMRG: The Warmup. 

The system, which we shall consider to be ID, starts with a small chain. 



This chain is split into a left and a right part, with two sites between them: 



Left Right 
This system shall be schematically represented as follows: 



In this figure two types of objects are to be seen: blocks (whether big, which contain many 
sites, or small, with a single site) and links between them, which we shall call hooks. Both shall 
be represented by certain operators. 

The DMRG growth loop works in this way: 

1. - Using the block and hook operators we compose an effective matrix which represents the 
superblock hamiltonian. 

2. - The ground state of this superblock hamiltonian is found. 

3. - With that state, density matrices are found for the left and right halves of the system: 



Left part Right part 

4. - We retain the m eigenvalues of the density matrix with highest eigenvalues on both sides. 
With them, appropriate truncation operators are built. 

5. - Using these truncation operators, the left and right parts are "renormalized" to yield new 
blocks and hooks: 



H--H 



6.- The hooks which arc dangling may not be tied up together, since they represent links to 
individual sites. We introduce two new sites between both blocks and we have closed the KG loop. 



After this process the system has increased the number of sites: N — > N + 2. 
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DMRG: The Sweep. 



The warmup by itself may yield rather satisfactory results: the diagonalization of the su- 
perblock hamiltonian provides at each step an approximation to the lowest energy states of the 
system. But, when combined with the sweeping process, the results of the DMRG may reach 
arbitrarily high precision^. 

The warmup provides a series of numerical values for the matrix elements of the block and 
hook operators which are coherent (in a sense which shall be now explained). All the operators 
which have been obtained through the warmup process must be stored. 

Let us call L(p) and R(p) the set of block and hook operators representing the p leftmost and 
rightmost sites respectively. In this way we may write the formal equation: 

L(p) + . + . + R(N-p-2) = n 

where Q represents the whole system, the symbols • represent individual sites and the addition 
symbol is suitably "extended" (see section 4.6). 

This very equation may be graphically represented: 

L{p) 

R{N - p - 2) 

Taking p — 2, we have 



L(2) + • + • + R(6) 



Left block 



Right block 



Since we have stored all the R(p) and L(p) during the warmup, it is possible to build up the 

superblock hamiltonian which represents the whole system. Once this has been done, the ground 
state of the superblock (target state) is found and the density matrix is built for one side. In our 
example, we take the left one: 



Lip) 



R(N-p-2) 

We obtain the most probable states (i.e.: the eigenstates of the density matrix with highest 
eigenvalues) and "renormalize" the pack L(p) + • so as to make L(p + 1) up. Afterwards, we search 
in the "warehouse" the pack R(N — p — 3), which "matches" the L(p + 1) we have just built to 
make up another superblock: 



L(p + 1) + . + . + R(N - p - 3) = Q 

In our concrete example, p takes now the value 3. Then 4, 5... 



This is true for QM. For interacting systems, the number m of conserved states must be high enough. 



3.1. The Density Matrix Renormalization Group. 
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Left Block 



Right Block 




L(3) + . + . + R(5) 
L(4) + • + • + R(4) 
L(5) + • + • + R(3) 



L(6) + . + . + R(2) 



This process goes on, making the left block grow and the right one shrink until the later is as 
small as possible. Once this point is achieved, we invert the sense of the sweep and the left block 
starts decreasing (and the R blocks arc the ones which get modified). 

Once the initial point is regained, it is said that a full sweep has been completed. Energies 
usually converge in a few sweeps up to machine precision. 

Some details were removed from the former description for the sake of clarity. 

• At each step only a given number m of states are conserved. But the packs L(p) and R(p) 
only contain the matrix elements of the needed operators between those states. 

• There is a minimum block size, which we shall term po, which depends on the number of 
conserved states m. In other words, po sites should contain exactly m states. 

• When the warmup procedure finishes, values for all L(p) and R(p) have not been obtained, 
but only for 

{L(po), L(po + 1), . . . , L(N/2 - 1)} {R(p), R(po + 1), . . . , R(N/2 - 1)} 

The reason is that the two last blocks of each series, along with two free sites, make up the whole 
system. Therefore, the sweep always starts with this configuration, making the left block grow 
(e.g.) until it reaches it maximum size L(N — Po — 2). 



The Key Features of DMRG. 

The DMRG is a high precision method. In the next section we shall describe its application to 
quantum mechanical problems and we shall discuss its performance. Before providing numerical 
results we would like to remark some of its basic features: 

• It is a variational method. There are no perturbative series: the system may be as strongly 
coupled as desired. Meanwhile not all implementations guarantee the lowering of the energy 
of all states at every step, it is sure that the accuracy of the approximation of the target state 
shall increase. 

• It is an implicit method: the states we are working with arc only known through their matrix 
elements of the operators making up the packs L(p) and R(p). Hopefully, the number of such 
values is much lower than the dimension of the Hilbert space. 
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3.2. DMRG Algoritms for Particles in a ID Potential. 

Now we shall describe a concrete non-trivial implementation of the DMRG algorithm: the 
solution to quantum mechanical problems in ID. This application was developed by M.A. Martin- 
Delgado, R.M. Noack and G. Sierra in 1999 [MSN 99]. 

Let us consider a similar problem to the one we found in section 2.2, i.e.: the problem of a 
particle in a discretized box with a potential. We write the matrix elements of the total hamiltonian: 

ifi = j 

Hij = < -l/h2 if = 1 

otherwise 

where N is the number of cells and h = Ax = 1/N, since the box is considered to be the interval 
[0, 1]. The values of Vi are "sampled" at each cell, e.g.: Vi = V((i — l)h,). 

• At each KG step there are some states \'^\) with i G [1 . . . m] for the left block and other 
for the right one. These states are not explicitly stored. The left block embraces p sites and 
the right one N — p — 2. The two central sites are represented by two delta states (concentrated 
at a single cell), denoted by "cZ" (center-left) and "cr" (center-right). 



> O O Q Q Q Q 



|5cr) 



Figure 3. States making up the superblock: \tpi'') correspond to the left block, \S^i) and \5cr) to the free central 
sites and jV'f^) corresponds to the right block. 



• Let US consider all possible decompositions of the system into a left part, two sites and a right 
part. Let us see how many packs we need: 

{L(m), L(m + 1) . . . L(N - m - 2)} {R(m), R(m + 1) . . . R(N - m - 2)} 

A block smaller than m sites makes no sense, since it would not even contain m independent 
states. 

• Each one of the packs, e.g. L(p) contains a block matrix and a hook vector. We shall denote 
the matrix with the same symbol (abusing of notation): L(p)tj and the vector ones by TL(p)i. 
In an explicit fashion: 

L(p),^ = (^^{-| H \^\>l) Mv)i ^ {^\\ H |6ci) [1] 
On the right side the definition is equivalent: 

R(p)tj ^ (ibf I H \^\>f) Mv)i ^ {^\>f\H \8,r) [2] 
These are the data we shall work with. 



3.2. DMRG Algoritms for Particles in a ID Potential. 
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Blocks Composition. 

One of the keys of the new renormahzation methods is its generalization of the blocks fusion 
idea. For the DMRG each block is represented by a "pack" , which consists of a series of matrix ele- 
ments of certain operators. Two packs may get fused and yield other pack which is a representation 
of the block including both of them. 

Let us suppose the system to be divided into a block with p sites, two free sites and another 
with q sites such that p + 2 + q = N . We write formally: 

L(p) + • + • + R(q) = 

In this paragraph we shall describe the process 

L(p) + .^L(p + l) 

which provides a representation (a pack) for the block containing the first p + 1 sites. 

Foremost, we shall write down a variational Ansatz which includes all the states of the su- 
perblock shown in figure 3. 

m m 

Id)) = al \^\>l) + a-' + a" |6er) + E l^f > [3] 

i=l j=l 

The variational parameters are the {<i| }]^-|^, {if a'^^ and a*^^: in total 2m + 2. Since the 
set of states is orthonormal, obtaining the superblock hamiltonian is equivalent to calculate the 
matrix elements of the hamiltonian between each pair. 

Using the matrix elements L(p)ij, R(q)ij, TL(p)i and TR(q)i we may make up the superblock 
hamiltonian like this: 



/ 



H 



Sb 



V 



L 








Tl 


Hcl,cl 


Hcl,CT 






^CT,cl 


Hcr,cr 


Tr 






4 


R 



[4] 



The diagonalization of such an effective superblock hamiltonian shall provide us with an esti- 
mate for the lowest energies of the total system. Its eigenvectors shall yield global states for the 
total system. 



5=1 5=1 

with b-^ denoting the j-th component of the i-th eigenvector of the superblock hamiltonian (j G 
[1...2m + 2], i G [l...m]). 



Now the properly called "renormahzation" procedure is executed, which we are denoting with 
the formal expression L(p) -I- • — >^ L(p -I- 1). 
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Let US focus on the left half of the system (left block + free site). We put zeroes on the 
components with j ranging from m + 2 upto 2vrL + 2 of matrix b and we obtain the states 

m 

5^bi|tl;!-) + b-+M5ci) 
i=i 

These do not make up an orthonormal set. A Gram-Schmidt process solves the problem and 
the correct states shall be expressed by 

m 
5 = 1 

With these new states, including one more site than the ones before, shall be used to build up 
the new pack L(p + 1). The block matrix elements shall be 



i',i' = l 



y=i 



i'=l 



This expression, apparently rather complicated, may be written as a basis change on the upper 
left part of the superblock hamiltonian: 





t; 


L(P) 







Bt 



[5] 



The renormalization of the hook is even more simple. Since it denotes the linking of the block 
to an external site, it is by definition 

\{V + l)i = 



H|6cr) 

but the whole block is only linked to the central-right site through the site which was formerly the 
central-left one (and is now part of the block). Thus: 



Mv + l)i = BT-+1 (Sell H |6cr) = BT^+iHp+i,p+2 



[6] 



To sum up, the process L(p) + • L(p + 1) is divided in these steps: 

• We take the left and right packs (equations [1] and [2]). 

• We make up the superblock hamiltonian (equation [4] ) . 

• It is diagonalized, and we retain the m first states. 

• We make zero the last m + 1 components of the vectors and re-orthonormalize them (Gram- 
Schmidt), to obtain a basis change matrix B. 

• We use that basis-change matrix to obtain, from the appropriate quadrant of the superblock 
matrix, the new L(p -|- 1) (equations [5] and [6]). 



3.2. DMRG Algoritms for Particles in a ID Potential. 
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Obviously, the scheme for the right side renormaUzation is fully analogous. 

N.B.: The density matrix is not necessary for the quantum mechanical problem. It is possible 
to introduce it, as it is done in [WHI 92], but the same result is obtained. Despite that, the method 
is so similar in spirit that it conserves the name. 

Warmup Cycle. 

In the first step {ab initio) the system is formally represented by 

n(2m + 2) = L(m) + • + • + R(m) 

The left block states, with m sites, may be (e.g.) the delta states concentrated on each of the 
m first sites. Mutatis mutandis, the same may be said about the right side. The "central left" and 
"central right" sites shall be respectively the sites m + 1 and N — m. 

Therefore, the packs for this first step are: 

L(m)i^ = (5i| H |6j) TL(m) = (6i| H |5^+i) 
R(m)ij = (6N+i_i| H |6N+i_j) TR(m) = (5N+i_i| H |6N-m) 

We renormalize it according to the formerly stated scheme both for the left and right sides: 

L(m) + »^L(m+l) R(m) + • ^ R(m + 1) 

by inserting two new sites between both blocks. The process is repeated until the total size of the 
system is reached: 

L(N/2 - 1) + • + • + R(N/2 - 1) = fi(N) 
At this very moment we state that the warmup is finished. All the packs 

{L(m) . . . L(N/2 - 1)} {R(m) . . . R(N/2 - 1)} 

have been initialized. Even though not all the L(p) and R(p) have been initialized, it is possible to 
start working with the values we have so far. Notice that, if the hamiltonian is symmetric under 
parity, it is possible to assume that, at every step, L(l) = R(i), and the problem gets simplified. 

N.B.: Obviously, there is no real link between the two sites "•+•" at every step of the warmup. 
It is necessary to include a fictitious one so as the algorithm works. When the hopping elements 
between neighbouring sites are equal all along the lattice, this is not a problem. 
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Sweeping Cycles. 

At a practical level, a sweep cycle starts with the system as the warmup leaves it. Although 
we might be tempted to start with the smallest possible left block (m sites) and the biggest right 
one (N — m — 2), the difficulty is obvious: only the blocks with sizes smaller than half the system 
are available. 

Therefore, the sweep consists of three parts: 



The left block is made grow against the right one till the extreme is reached (from N/2 — 1 to 
N -m-2). 



The right block now grows from its minimum size m to its maximum N — m — 2. 



The left block grows again from m to N/2 — 1. 



Once the third phase has been finished, a sweeping cycle has been completed. All the blocks 
which include less than N/2 sites have been refreshed at least once from their creation and the rest 
of them have been created in the process. 

The DMRG algorithm repeats the sweeping cycle until convergence is attained. If the values 
for the superblock hamiltonian eigenvalues are not modified (within some precision) in a whole 
cycle, the problem is assumed to be solved. 



Numerical Results. 

In [MSN 99] the behaviour of the preceding algorithm is analyzed for a variety of potentials, 
and results are obtained whose errors are smaller than one part in 10^° for all of them: harmonic 
oscillator, anharmonic oscillator and double well. 

Maybe the most interesting results were obtained with the double well potential, where in 
perturbation theory it is rather difficult to obtain the modification of the spectrum due to the tunnel 
effect between them. This modification is analytically studied perturbatively using the instantonic 
approximation. The results which are obtained are fully consistent with exact diagonalization 
whenever possible and with the instantonic approximation when it is not. 

Since these numerical results were not obtained originally for this work and are not relevant 
for its understanding, we refer the reader again to [MSN 99] for its careful analysis. 

It is important to remark that the computation is rather fast because the process is implicit, 
i.e.: the components of the complete wave- vectors are never explicitly used, and they are not even 
stored. The matrix elements of specific operators between those states is enough to carry out the 
RG cycle. 

This feature of being implicit is quite desirable for an RG algorithm based on the block idea, 
since the number of operations to take an RG step does not scale with the system size. Unfortu- 
nately, as we shall see in the rest of this work, this is not always possible. 



3.3. DMRG for Trees. 
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Part II. DMRG, Trees and Dendrimers. 



3.3. DMRG for Trees. 

The implicit DMRG algorithm which was described in the former section was strongly depen- 
dent on the left-right distinction. Is it possible to find an algorithm with the same features on an 
space where this distinction does not exist? 

A tree is a graph which contains no loops. In a certain sense, it keeps the unidimensional 
character which is needed^ for DMRG and, at the same time, it prepares the terrain to tackle 
afterwards multidimensional problems. In [MRS oo] M.A. Martm-Delgado, G. Sierra and the author 
of this thesis extended the quantum-mechanical DMRG algorithm so as it might work on a tree-like 
graph without losing its implicit character. This technique found a practical application in the 
study of excitons on dendrimers, which is developed in the next section. 

Climbing the Trees. 

In the former chapters we have defined the hamiltonian operators on discretized spaces based 
on a graph structure. In this section graphs acquire a crucial importance. 

Let ^ be a graph, described by a set of sites S and an associated neighbourhood structure: for 
each i £ S there is a set N(i) C S such that i £ N(j) ■^=^ j G ^{^)- usual, a path shall be an 
ordered set of sites {ii . . . in} such that i^+i G N(ii^) (i.e.: each site is connected to the previous 
one). The number of elements of N(i) is called d(i), the degree of the sites or coordination index 
in physical terminology (see appendix A and [BOL98] for more details). 

On a connected graph, the removal of a site i splits the set S into a certain number k(i) of 
disconnected components: S(i)i . . . S(i)]^(|), and it must always be true that k(i) < d(i). If the 
inequality saturates for all the sites, the graph is called a tree (see figure 4). 




Figure 4 A tree: the removal of any site leaves the graph split in disconnected components such that each 
neighbour of the site remains in a different one. 

It is important to consider the previous definition in more intuitive terms (however, this ex- 
planation does not make up a proof, and we refer the interested reader to [BOL98]). 



There are approaches to tackle 2D systems with DMRG (see next chapter), but the error range is much 
higher and the speed much lower than for the ID case. 
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Let us consider a connected graph in which, when removing a site p which contains three 
neighbours, there appear three disconnected components. It must be true that each of the com- 
ponents contains one of the neighbours. Let us call the three components S(p)i, S(p)2 and S(p)3. 
If the site p did not exist, it would be impossible to trace a path from a site in S(p)i to a site of 
S(p)2- Thus, the path between those sites passes necessarily through p. 

If that is true for any p in the graph, then it is not difficult to conclude that for any pair of 
points in Q there is just one possible path and, therefore, there are no closed paths (which imply 
the presence of two different paths between two points). This is the classical definition of a tree. 

As a matter of fact, there are various alternative definitions of a tree, and it may be useful to 
know them: 

• A connected graph for which the removal of any site leaves the graph split into as many isolated 
components as neighbours it has. 

• Connected graph for which the removal of any site leaves its neighbours disconnected among 
themselves. 

• Graph in which for each pair of sites there is a single path which connects them without 
repeating sites. 

• Connected graph without loops (i.e.: there are no closed loops which do not repeat sites). 

It must be clearly remarked that "legitimate" multidimensional spaces may not be represented 
as tree-like graphs. In some way, trees present extended unidimensionality. On the other hand, 
a unidimensional space with periodic boundary conditions is not a tree. The final thesis of this 
chapter shall be that "The implicit DMRG technique finds its natural place in the trees" . 

DMRG ON Trees. The Idea. 

Let us see how the DMRG algorithm for ID chains is adapted to tackle problems on trees. 

It is fundamental to observe the system from the "point of view" of a particular site. Let p be 
that site. For it, the universe is reduced to a certain number of blocks (d(p) to be exact) which 
are disconnected among themselves, as it is shown in figure 5. 




Figure 5. The world as it is seen from site p. 

If we are endowed with appropriate "packs" (block + hook) for each of these blocks, we may 
build a superblock hamiltonian which represents the whole system. This hamiltonian would be 
diagonalized and, with the help of its m lowest energy states, we would carry out the block fusion 
process. For example, blocks S(p)i and S(p)2 along with the site p itself may form a block 

S(p)i + S(p)2 + »p 

Graphically this appears in figure 6. 



3.3. DMRG for Trees. 
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Figure 6. Blocks S{p)-i and S{p)2 along with site p get fused to make up block B. 

Let us call q the neighbour of p which belongs to the block S(p)3. Block B is a part of the 
system seen from q. Thus, there must exist an index y such that 

B = S(q)y 

Thus, the operation is really 

S(p)i + S(p)2 + »p ^S(q)y 
Now we may consider site q to be free. 



A Simple Illustrative Computation. 

Let us consider the graph which is represented in figure 7. 




Figure 7. A sample tree, on which we shall apply the DMRG algorithm. 



Let site 2 be our free site for the first step. The system, as it is seen from that site, consists of 
two blocks (besides itself). The first contains only site S(2)i = {1}, while the second contains the 
sites 5(2)2 = {3,4,5,6,7}. Figure 8 shows this graphically. 




Block S(2) 
Sites: 3.4 




Figure 8. Schematic representation of the system viewed from site 2. 



Let us suppose that only the ground state of the system is to be found and we associate a 
single normalized state to each block. For the first one we have the delta state |i|)(S(2)i)) = |6i). 
On the other hand, the state for the second block, |tl;(S(2)2)), encapsulates 5 degrees of freedom 
(as many as sites) in a single one. Thus, we may estimate variationally the ground state of the 
system using as an Ansatz 

= ai |tKS(2)i)) + a2 |t|)(S(2)2)) + I62) 
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The superblock matrix, which "represents" the full system, is a 3 x 3 matrix which takes this 
form the point of view of site 2: 



H 



Sb 



H[S(2)i] T[S(2)i] 

H[S(2)2] T[S(2)2] 
Jt[s(2)i] Tt[S(2)2] H22 



Since the three states are normalized, the hamiltonian matrix elements between them yield 
the elements of H^"^. 



H [S(2)i] = (tl;(S(2)0| H |xl)(S(2)i)) H [8(2)2] = (tKS(2)2)| H ^{2)2)) 
T [S(2)i] = (iKS(2)i)| H I62) T [S(2)2] = {MS{2)2)\ H I62) 

Once we have obtained the ground state we may carry out the blocks fusion procedure 



S(2)i + .^B 

Graphically, this is shown in figure 9. 



{1} + .^{1,2} 




Figure 9. Block renormalization: the old block S{2)i plus the free site 2 get renormalized to the new block B. 

If the following free site is 3, then the block B corresponds to one of its constituents, let us 
suppose it is S(3)i. Let us see how to calculate its block matrix element and its hook. We eliminate 
from state |4)) any mention of the block 8(2)2 and renormalize: 

|il;(S(3)i)) = Ai |il)(S(2)i))+Ae|62) 
where Ai = N • ai and Ac = N • Qc, where N is an appropriate normalization constant. Now, 



H[S(3)i] = (•il;(S(3)i)|H|tKS(3)i)) = |Aip H [S(2)i] + |Ac|' H22 + 2Re (AiAcT [S(2)i]) 

and now only the "hook" is missing 

T[S(3)i] = (iKS(3)i)|H|63) = AcH23 

Once we have completed the process, site 3 may be set free, which will have blocks S(3)i = 
{1,2}, 5(3)2 = {4, 5} and S(3)3 = {6, 7}. The pack for the first of those blocks has been determined 
at the previous step, and the other two packs must be taken either from a former sweep (RG cycle) 
or from the warmup. 



3.3. DMRG for Trees. 
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S(qX 



S(P)> 



Figure 10. Graphical description of blocks S{p)x and S{q)y. 



Blocks Composition. 

The major difference between ID DMRG and its adaptation to trees is the blocks fusion 
technique. We now give the general rule. 

Let us suppose that, at a given RG step, the free site is p, while the following one shall be q. 
Obviously, it must be true that q G N(p). 

Let S(p)x be the block from p which contains its neighbour q. And let S{q)y be the block 
from q which contains p, as it is reflected in figure 10. 

The objective, therefore, is to carry out the formal prescription 

J2 S(p)i + .p ^ S(q)y 
ie[i...d(p)]-x 

Le.: the "addition" of all the blocks corresponding to site p, but the one labelled x, along with 
the site p itself to "refresh" the block from site q which contains p. 

Each block possesses a series of m states which constitute an orthonormal set^. Let us call 

the l-th state (I G [1 . . . m]) of the k-th block when sphtting the graph according to site p. With 
these data we shall constitute the variational Ansatz (notice how this expression generalizes [3] ) . 



d(p) m 

I*) = X! 0""^""^^+' K(S(p)k)'> + C"^'^(T')+1 
lc=l i=l 



[7] 



The data pack is now formed by the set 



H[s(pKkj ^ (iKS(pKr|H|iKS(pK)'> 

T[S(pK]i = (iKS(pKr|H|6p) 

With these elements the superblock hamiltonian is easy to construct, and has a very charac- 
teristic structure 



This number might in principle vary from one block to another, although we have considered it to be uniform 
for the sake of simplicity. 
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H 



Sb 



V 



H[S(p)i] 




T[S(p)i] 




H[S(p)2] 




T[S(p)2] 










TMS(p)i] 


[S(P)2] 







[8] 



This superblock hamiltonian needs to be diagonalized. The eigenvalues are the actual estimate 
to the energies of the full system. The m lowest eigenstates, making up the matrix Bj- with 
i G [1 . . . m] and j G [1 . . . md(p) + 1] undergo the following process. 

• Making the components corresponding to block S(p)x zero, i.e.: from (x — l)ra+ 1 to xm. 

• Re-orthonormalization of the vectors, through a Gram-Schmidt process. 

We shall call the same matrix after having finished both processes. Now we may assert that 
the new states which correspond to block S(q)y are 



d(p) 



m(k-l)-|-i 



|il;(S(pKr) + B"^'i(T')+i|5^) 



k=l i=l 



The fact that the summation extends also to k = x must not bother us, since the corresponding 
values of B are null. 

With the new block states ready, we proceed to renormalize the block matrix and the hook 
vector of S(q)y : 

md(p)-|-l 

H[S(q)Jtj = (;I.(S(q)y)^|H|4)(S(q)y)')= J] BIbH\ 

lc,l=l 

The hook requires the matrix element between sites p and q. The probability of presence in 
site p is given just by B^^^^^^^. Thus, 



T[S(q)y]i = ^]nd{p)+l^V^ 



3.3. DMRG for Trees. 
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The Warmup. 

Once we have analyzed the tools, we tackle the concrete algorithm. The first step we must 
take is the obtention of a sufiicient number of packs so as to start the process. 

The only condition on the initial packs is their coherence. We mean that it is necessary that 
there exist real states whose matrix elements arc given by their values. If the initial packs are 
coherent, the distance of the initial states to the true ones does not matter: convergence is ensured 
(faster or slower). 

It is impossible to give an universal algorithm for the warmup without lowering the efficience 
in particular cases. The last section of this chapter deals in detail with a process for the analysis 
of dendrimeric molecules. For the moment we state some general rules. 

If the tree is not infinite, there must exist a set of sites which only have one neighbour. Let 
Di = {p G ^|d(p) = 1}. Once the computation of those sites is finished, they are removed from the 
graph. Let Gi be the resulting graph after that elimination. Now we define D2 = {p G ^i|d(p) = 
1}. This process yields a sequence of sets Di^ which gives the different onion skins of the graph. 
The number of steps to take until there are no more sites in the tree shall be known as its depth. 

The warmup algorithm must proceed inwards, i.e.: shall start by grouping into blocks the sites 
belonging to the lower order onion skins, which shall absorb the inner sites in an iterative way until 
they reach the center. At that moment, the sweeping cycles may start. 

An example may be useful in this case. Figure 11 shows the first step on a tree. The sites of 
the first and second onion skins are put into blocks. If these blocks contain more than m states, it 
is necessary to truncate. Otherwise, all the states are conserved. 




Figure 11. Warmup of a part of a tree. The arrow on site number 1 shows that the graph goes on in that 
direction. The sites bagged in the left side belong to the first and second onion skins. On the right, some blocks 
appear which have substituted them, given in parenthesis the number of states they contain. If this number is 
< m, it remains like that. If it is greater, it is truncated to m. 



The warmup is complete when all the sites belong, at least, to one block. At this moment the 
sweeping may start, and at its first step the free site must be one of the ones which belong to the 
innermost onion skin. 
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The Sweeping Cycle. 

The strategy starts by the estabhshment of a sequence of free sites along the tree which we 
shah call the path. Generally speaking, the path may visit each site more than once. 

This sequence does not need to cover the whole graph. It is not necessary to enter into a block 
when the number of sites which it contains is < m (i.e.: there is no truncation). For example, if 
m = 3, it would not be required to enter at figure 11 into any of the blocks which are marked. 

The election of the path, and that of the warmup, must be built taking into account the 
concrete problem, so as to take profit of the available symmetries. The case of the dendrimers shall 
provide a clear example. 



3.4. Physics and Chemistry of Dendrimers. 

A practical application of the DMRG technique applied to trees is the study of the excitonic 
spectrum of the polymeric molecules called dendrimers. The main reference for that analysis is the 
work of M.A. Martm-Delgado, G. Sierra and the present author [MRS oo]. 

This section is dedicated only to provide a good physical and chemical basis to the dendrimers' 
analysis. A great amount of basic information may be obtained from [VOG 98] and the references 
therein. Along this section a part of the classical works on the subject are cited, along with some 
of the recent developments which are relevant for our work. 

Dendrimers. 

Dendrimers^ are highly branched polymers whose size in some cases exceed 10 nanometers. 
They are synthesized molecules with a rather low polydispersity (i.e.: the molecular masses dis- 
tribution is highly peaked around a value) which have quite peculiar architectural features. These 
are: 

• They have a "core" , called the focal point. 

• Out of this core stem a number of branches, called wedges. 

• After the addition of certain number of monomers, equal at all the branches, these are forced 
to branch again. A generation has been completed. 

• Successive generations are formed by adding a fixed number of monomers to each branch until 
they are forced to branch again (see figure 12). 

The number of generations is not unlimited. The molecules fold up in three-dimensional space 
(they are not flat), but steric repulsion (exclusion volume effects) inhibits high generation numbers. 
The exact saturation generation depends highly on the nature of the monomers and the molecular 
architecture. 

So as to simplify the analysis it is a common practice to ignore the chemical properties of the 
monomers, which are now considered to be merely "sites" or "nodes". We might say that we embed 
the molecule into a lattice or a graph for our purposes. 



The word dendrimer comes from Greek 5kv5pov, which means tree and fiepo<;, which means part. The term 
was coined after its resemblance to dendrites. 



3.4. Physics and Chemistry of Dendrimers. 
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Figure 12. A dendrimer with 3 wedges, 4 generations and connectivity 3 at the branching points. 

Chemistry of Dendrimers. 

The synthesis of dendrimeric molecules was under very active search by the end of the 70's, but 
it became realistic when Fritz Vogtle's concept of molecular chain reaction [B w V 78] was developed. 

Vogtle's idea, known nowadays as the divergent method of synthesis, works by adding iter- 
atively shells from the central core. Chemists may control the process through the addition and 
removal of molecules which allow or inhibit polymerization and branching. Through this process a 
high control on the molecular architecture is gained, but is rather difficult to carry out in practice 
[TNG 90]. 

A new method was developed by Hawker and Frechet [HF 90] , which is called the convergent 
method, which makes dendrimers grow from the surface inwards. The small pieces which are 
formed are afterwards "self-assembled" until the molecule is complete. The great advantage of 
this second technique is its mechanizability. 

Dendrimeric molecules may employ a great variety of monomers. Typical elements are amino- 
acids, polyamidoamines, DNA bases or glucids. Recently also organo-silicic and other organo- 
inorganic hybrids have been used. 

Of special interest to our work are dendrimers whose monomers are phenyl rings, which are 
linked by acetylene molecules. This link is usually termed diphenylacetylene. The molecules are 
named as Dn, where n is the number of phenyl rings (e.g.: D4, D127...). 

Which are the good properties of dendrimeric molecules which made them so actively pursued 
by the scientific community? 

The geometrical configuration of the molecule makes it have many "free ends", which may 
be functional. Usual polymers with a linear chain have only two such terminal monomers, and 
traditional branched polymers have a number of them which is a small fraction of the total. In 
dendrimers, on the other hand, they suppose a non-negligeable fraction of the total (even beyond 
one half). 

Its self-similar structure makes it have a large surface to volume ratio [GEN 79]. Concretely, 
they have been considered for their usage as "recipients" , where guest molecules may travel through 
hostile media (e.g.: drugs deliverers). These hosting systems are based on the compact packing of 
monomers on their surface if the saturation generation has been achieved. It has been experimen- 
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tally tested [JBBM 94] that trapped molecules have a neghgeable trend to diffuse out of the box, 
but the interaction between recipient and guest is still under discussion. 

Another application which is also based upon their fractal structure is the possibility of building 
artificial molecular antennas for the infra-red range [KOP 97]. The wavelengths which an antenna 
may deflect depend on the lengths of the "corridors" that may be travelled by the electrons undis- 
turbed. Dcndrimcric molecules may be built so as their corridors range a wide spectrum of lengths. 
We have focused on this application. 

Let us cite at last some applications which have either been already carried out or are being ac- 
tively researched: gene therapy, nanotechnology (as building blocks to construct bigger molecules), 
magnetic resonance imaging (as contrast agents), molecule recognition, antiviral treatments, en- 
ergy converters (of light into electric currents) [VOG 98]. A host of applications which range from 
materials science to medicine. 

EXCITONS IN DeNDRIMERS. 

As it was previously mentioned, one of the applications of dendrimers which may attract the 
attention of physicists is the possibility of building antennas for harvesting infra-red light. It is 
necessary, therefore, to model the behaviour of the electrons. 

The Frenkel exciton model is a development of the tight binding model (TBM) for an electron 
in a lattice of similar atoms or molecules, which we shall merely term sites. The approximation 
takes orbitals at each site as its starting point. 

Let us suppose that there are only two orbitals per site which contribute to the physics of the 
problem. We shall call |0i) the orbital of site i with smaller energy and |li) the excited one at the 
same site. 

Let us introduce a vacuum state per site |vaci), meaning that there is no electron at all at the 
site. A couple of creation and annihilation operators is suitable: 

c{o K^^i.) = |0i) ct^ |vaci) = |li) 
As it was expected, Coi and Cn shall denote the corresponding annihilation operators. Let us 
define now excitation operators: b| = c||Coi and bi = cJ^Cii. These operators fulfill 

b||Ot) = |li) b||lt)=0 
bi|Oi)=0 bi|lt) = |Ot) 

Le.: the excitation operators act as creation and annihilation of a new "quasi-particle" : the 
exciton. 

The most simple hamiltonian which may be written for the dynamics of these excitons is 

ies (ij> 

where S represents the set of sites of the molecule, Et are the energy gaps at each site (they may be 
different) and play the role of chemical potentials. The sum extended over the (i, j) runs, as it is 
usual, over the graph links. The origin of these terms is the overlapping of the orbitals at different 
sites, which may make the electron "jump" from one site to another. 

An electrostatic interaction might be added between different electrons to obtain a much more 
complex interacting hamiltonian such as, for example, the Pariser-Parr-Pople hamiltoinian (which 
has been analyzed with DMRG [RT99] and similar techniques [MSP J 99] for ID chains). 



3.5. DMRG Algorithm for the Dynamics of Excitons in Dendrimers. 
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Since the model formerly introduced [9] is free, it may be converted into a problem in quantum 
mechanics just by writing down a wavefunction 

\A>) = Y,Ci\h) 

ies 

And, therefore, in the one-particle sector, the hamiltonian gets reduced to 

Hiex = + OiA (6j I + H.c.) [10] 

tes (i,j> 

Schrodinger's equation is 

and the problem has been reduced to another one analyzable with DMRG techniques for quantum 
mechanics. The eigenvalues of Hiex, anyway, must be taken carefully. The lowest eigenvalue does 
not mean the energy of the ground state, but the minimum energy which is needed to create an 
exciton. In the terminology of this research field, it is called the optical absorption edge. 



3.5. DMRG Algorithm for the Dynamics of Excitons in Dendrimers. 

The physical problem which we shall study is the prediction of the optical absorption edge for 
the series of dendrimers of the diphenylacetylene family. There are experimental results, obtained 
by Kopelman et al. [KOP 97] to which we shall compare our results and those of other authors. 

The theoretical model we shall employ was developed by Harigaya [HAR98] [HAR99]. It is 
a modification of the excitonic model given in the previous section, with disorder added at the 
nondiagonal terms. 

Our target shall be, therefore, the analysis of the hamiltonian [10] through DMRG. A set of 
realizations for the set of values {Jij} shall be considered and the optical absorption edge shall be 
obtained statistically. 

Even though it is well known, it is appropriate to remark that quantum-mechanical models 
with disorder are equivalent to many-body interacting models [ITZD89]. The process consists 
(briefly summarized) in the conversion of the summation on realizations into a path integral for a 
new fermionic field. 

Compact and Extended Dendrimers. Combinatorial Properties. 

A Cayley tree (also known as Bethe lattice) is a tree with a uniform degree d(i) = c. It is 
impossible to make a finite model of such a tree, so we may define a finite Cayley tree as that 
which has only two kind of sites: terminal with d(i) = 1 and normal with d(i) = c. 

A dendrimer whose sites are linked describing a finite Cayley tree is called a compact dendrimer. 
The connectivity is usually c = 3 (see figure 13), since it is the smallest non-trivial value. We shall 
only consider dendrimers with complete generations. 

The compact dendrimer of figure 13 has 4 generations and shall be denoted by Dq[4). The 
number of sites of such a graph (translating: the number of monomers of such a dendrimer) with 
L generations is given by 
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Figure 13. A compact dendrimer, i.e.: a finite Cayley tree. A further generation would force some sites to 
overlap (or the molecule to fold up in 3D), which exemplifies steric repulsion. 



They are also termed with the notation Dn, where u is the number N(L). Thus, the first molecules 
of the compact series are D4, DIO, D22, D46,... 

An extended dendrimer is something more complex to describe. Generations have different 
lengths. The branches which are close to the core are longer than the ones which are near the 
surface. If L is the number of generations and n is the generation index (counting outwards), the 
length of the branches is L — u. A "termination shell" is added, which would correspond to u = L. 

As an illustrative example, figure 12 may be good, which has L 4 generations. The innermost 
one, u = 1, has length L — u = 3, and at each branching it decays by one site. 

The number of sites of an extended dendrimer may be analytically found. The computation 
is as follows: 

• A central site. 

• The three first branches: in total 3(L — 1) sites. 

• The second level branches: in total 3 • 2 • (L — 2) sites. 

• The n-th shell gives 3 • 2^~^ • (L — n). The shells are considered from 1 upto n — 1, since there 
is no sense in adding a null shell. 

• The terminating shell: 3 • 2^~^. 

Summing up: 



L 




n=l 



L 



Ne(L) = 1 + 3 J]](L-n)2^~^ + 3-2 



L-1 



u=l 



The final result is (see [GKP 89], where techniques of discrete analysis are discussed): 



Ne(L) = 1 + 3- (3 • 2'--^ - L - 1) 



[11] 
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The advantage of this concrete form is that it makes exphcit the contribution of each wedge 
plus the central site. Thus, the first dendrimers of the extended series are D4, DIO, D25, D58, 
D127,... Notice that the first two elements coincide in both series. 

We remark also the fact that the function N(L) presents an exponential growth in both cases. 
For relatively high values of L, the exact diagonalization of the hamiltonian is impracticable. 

Dendrimers and Computation. 

So as to make the numerical computations easily reproducible, we now give a series of useful 
technical details. 

The first practical problem is to write down a program which implements computationally the 
neighbourhood structure. Once this has been computed, the calculation of the hamiltonian matrix 
is straightforward. 

Although the programs are written in C++, a typical structure of the LOGO language was quite 
useful [AdS8i]. It is the turtle: an object which is able to advance and rotate and, as it moves, 
trace a drawing. Our turtles are also able to reproduce and, thus, with the collaboration of the 
whole family, draw the graph of the hamiltonian. 

We give now in pseudocode the instructions which we have used, both to draw compact and 
extended dendrimers. The routine of pseudocode 1, called Launch, only receives a parameter, 
which is a turtle. This has two inner variables: the generation number and the direction it is 
pointing at. 

The program invokes itself using as parameters each of its child turtles. It is, therefore, 
recursive, and the finalization condition is to have reached the last generation. 



Advance . 








Make new site. 








level ^ level — 1 . 








if (level > 0), 








Make new turtle with 


level 


level 


here . 


Rotate clockwise a. 








Launch . 








Make new turtle with 


level 


level 


here . 


Rotate anticlockwise 


a. 






Launch . 









Pseudocode 1. Instructions for a turtle to create a compact dendrimer. 



The complete program just creates three turtles (with angles 0, 2n/3 and An/3) and launch 
them. The PostScript graphics which are shown in this thesis were obtained with this program. 
Furthermore, the differences with the necessary code to create the neighborhood structure are 
minimum. 

In the case of the extended dendrimers a small modification is necessary: the instruction 
Advance must be carried out not just once but as many times as the length of the branch we are 
at. Pseudocode 2 gives the sequence of orders to substitute Advance. 
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For i = 1 to maxJevel — level 
Advance 
Make new site 

Pseudocode 2. Modification of the former program to create extended dendrimers. 

The instruction Make new site adds a vertex to the graph, assigns it the first available index 
and links it to the previous site. This automatic manner to generate the numbering of the sites of 
a dendrimer keeps an interesting relation to number systems. 

The key concept is the chain of decisions which must be taken so as to localize the site in the 
dendrimer. We shall consider the compact case for simplicity. 

At the central site four decisions are possible: one for each of the three branches plus the one 
corresponding to finish there. At each of the other sites there are only three possible decisions: 
"left" , "right" or, again, terminate. Thus, the specification of any site may be given by a word. 
The first letter belongs to the set {0,1,2,7}, and the following ones to {0,1,7}. Of course, the 
symbol 7 must be the last one in each word. Valid examples are 2100107, 107 or, simply, 7. 

Let us define now a bush of level n as the set of sites which would be obtained by launching 
a turtle at that level. The number of sites of a bush of level n is 2'^"'^+^ — 1. With this rule, the 
conversion of a word into a site index is quite simple: each numerical symbol contributes its value 
multiplied by the size of the bush of that level. The symbol for termination contributes a single 1. 

In a single formula: 

k 

{siS2...Sk7}^ 5^ Si •21-^+1 + 1 

i=l 

It is easy to convince oneself that this formula expresses a particular number system. 

The determination of a minimum path along the dendrimer is also an interesting problem. Of 
course, it is impossible to find a hamiltonian path (i.e.: one which runs over all the sites without 
repeating any one). The following algorithm returns a list with the order of the minimum length 
path which spans the whole tree. 

We shall start at the central site. The rule is "each time an option appears, choose the site of 
smallest index among the ones which have not been visited yet. If all the sites have already been 
visited, return the way you came from." 

Given any tree, the number of steps of the minimum path is the double of the number of links 
of that tree. Denoting this number by |P| we have 

N 

|P| = Ed(i) 

i=l 

This sum may be determined analytically for compact and extended dendrimers. In the first 
case we only have two kinds of sites regarding the number of neighbours: "inner" , which have 3 
neighbours, and "terminals", which have only 1. Therefore, the length of the path would be: 

|p|= E 1 + 

Terminals Inner 

The number of inner sites is nothing but the number of sites of the dendrimer with one 
generation less. Thus, 
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|Pc| = 3Nc(L - 1) + (Nc(L) - Nc(L - 1)) = 3 • 2^+^ - 6 

(let us remind that Nc(l-) is the number of sites of a compact dendrimer of L generations). In 
the case of the extended dendrimer the computation is shghtly longer. Sites may be classified 
according to the following criterion: 

• "Branching" sites, with 3 neighbours. 

• "Bridge" sites, with 2 neighbours. 

• "Terminal" sites, with 1 neighbour. 

The number of branching sites is equal to the size of the compact dendrimer with the same 
number of levels: 3 • 2^~^ — 2, and the number of terminal sites is the same as in the compact 
case: 3 • 2^~^. The number of "bridge" sites may be obtained by substracting from [11] or either 
by direct computation. Choosing the last option we have the summation 

L-2 

# "Bridge" sites = 3 • ^(L - i - 1) • 2^"^ = 3 • (2'--^ - L) 

i=l 

(notice that, although the formula was designed for L > 2, the final result takes into account 
correctly that for L = 1 and L = 2 there are no bridge sites). Putting all pieces together, the final 
result is 

|Pe| =9-2^-6-(L + 1) 

The Laplacian on a Dendrimer. 

Regarding only the connectivity properties, a dendrimer has a very large symmetry group. 
This group is Z3 (S> {^2)^, where n is the number of branching points of the molecule (excluding 
the focal point). It is possible to classify the eigenstates of the laplacian on a dendrimer, either 
extended or compact, according to the group representations. 

The ground state is, logically, uniform all over the graph and corresponds to the trivial rep- 
resentation of the group. Figures 14A and 14B show other eigenfunctions of the laplacian with 
free b.c. for an extended dendrimer with 6 generations which correspond to higher dimensional 
representations. 




Figure 14a. The two first excited states of the laplacian with free b.c. on an extended dendrimer with 6 
generations. Both of them have the same energy and correspond to a bidimensional representation of the molecular 
symmetry group. 
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Figure 14b. On the left, the seventh eigenstate, with sixfold degeneration. On the right, the eigenstate number 
166, which corresponds to an unidimensional representation. 

Harigaya's Model. 

The model proposed by K. Harigaya for the computation of the optical absorption edge of 
extended dendrimers of the diphenylacetylene family [HAR 98], [HAR 99] has the following features: 

• Each site (molecule) is represented by a single pair of molecular orbitals. 

• It is a quantum-mechanical model of a single exciton. 

• It has uniform chemical potential E. 

• The hopping terms between sites are random, drawn from a gaussian distribution with mean 
zero and standard deviation J. 

Analytically, the hamiltonian is given by 

N 

H = 5] E |5t) (6,1 + (Si I + H.C.) 

and the empirical constants are given by E = 37.200 cm~-^ y J = 3.552 cm~-^. 
Application of the DMRG Algorithm. 

The general algorithm exposed at section 3.3 is directly applicable to this problem. The 
warmup establishes the division into onion skins, which are translated in this case to the division 
into generations of the compact dendrimers. In the case of the extended dendrimers, generations 
are subdivided into different shells because of the corridor sites. Basically, there are only two 
schemes of blocks fusion: 

Branching site: Block + Block + • ^ Block 
Bridge site: Block + • ^ Block 

For our concrete application we chose m = 3, so our last two generations may be included in 
the first warmup step. Figure 15 shows the typical scheme. 

As the block goes "creeping" along a corridor, it "eats" the sites as DMRG does in ID. When 
two blocks which are going upwards arrive at the same time to a branching point, they are forced 
to get fused and keep ascending together. Figure 16 shows these two processes. 
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Figure 15. The last two generations are encapsulated into blocks. 




Figure 16. On the left side, a branching site is going to be fused with its two children blocks. On the right side, 
a "bridge" site (2) is going to be swallowed by the ascending block. 

Once the central site is the only one which does not belong to any block, the warmup may be 
said to be finished. The sweeping proceeds as usual. Each block is represented by m = 3 states, so 
the dimension of the hamiltonian matrix shall be 3ra + 1 = 10 at branching sites, and 2m +1 = 7 
at bridge sites. 

Numerical and Experimental Results. 

The target of our calculations is the optical absorption edge for extended dendrimers, as it was 
formerly stated. The procedure may also be used to study compact dendrimers, with the results 
which shall follow. 

The following table compares the experimental results with two theoretical predictions of 
different authors. 



Generations 


Experimental^ 


Theor. Kopelman ^ 


Theor. Harigaya ^ 


1 


31500 ± 100 




31511 


2 


31300 ± 100 


31250 


28933 


3 


27700 ± 200 


27244 


27750 


4 


25000 ± 300 


25800 


26000 


5 


25000 ± 600 


25300 


25022 



Table 1. Optical absorption edges measured in cm . (1) Experimental measurements of Kopelman et al. [KOP 
97]. (2) Theoretical predictions, also by Kopelman et al. [KOP 97]. (3) Theoretical results of Harigaya [HAR99]. 

The predictions of Kopelman are based on a thoroughly different theory, called linear cluster 
theory. This theory is based on the analysis of "corridors" and, according to it, the optical 
absorption edge is dominated by the longest length of a linear chain of monomers. The energies of 
the excitons travelling through them are calculated through Hoshen-Kopelman's theory [HK 77]. 
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Harigaya's model, which has already been exposed, is a model for a single particle. One of the 
goals of our study is to know whether it would be appropriate to undertake a more sophisticated 
study based on a many-body theory, such as the Pariser-Parr-Pople hamiltonian. The DMRG 
does not provide a new theory about excitons in dendrimers. It just makes it easier to obtain the 
numerical predictions of a given theory, if this theory fulfills some pre-requisites (as Harigaya's 
model does). 



Harigaya (theor.) o 
Kopelman (exp.) 

DMRG □ 



DMRG 




12 14 



Figure 17. Optical absorption edge experimen- 
tally obtained by Kopelman (as a function of the 
generation), along with the theoretical results of 
Harigaya and those obtained with DMRG. 



Figure 18. Optical absorption edge with its er- 
rorbars for Harigaya's model up to 13 generations 
obtained with DMRG. The fit to a power law is 
shown. 



Figure 16 shows the optical absorption edge against the generation number according to exper- 
iments, Harigaya's calculations (exact diagonalization) and the results obtained with DMRG. As 
it was to be expected, the last two coincide. The DMRG data were obtained with 10.000 samples. 

Our calculations may go beyond, and figure 18 shows the optical absorption edge up to 13 
generations along with their error bars. In the last two cases only 1.000 samples were drawn. The 
data fit reasonably well to a power law 

Ec(n) ^ 32000- Cen°'= 

with Ce ~ (2020 ± 200)cm-i and oie ~ 0.71 ± 0.04. Therefore, it may be asserted that the system 
does not have a "gap" in the thermodynamic limit. 

Since the DMRG calculations were carried out with m = 3, the next two states of the exciton 
were also calculated. Figure 19 shows them. Figure 20 makes us watch the eff'ect of the addition 
of wedges, which is a usual technical process, from w = 2 to w = 6. 

When we consider compact dendrimers, the results diverge greatly from the experiments. 
According to the data of Kopelman et al. [KOP 97] for the compact dendrimers series, the optical 
absorption edge hardly decreases when the number of generations increases. Concretely, Ec(rL) ~ 
31500ib200 for n € [1 ... 5]. I.e.: it is probable that the system has a finite gap in the limit n — > 00. 

Notwithstanding, when analyzing Harigaya's model, the fit to a power law of the results is also 
rather accurate: 

Ec(n) « Eo - Ccn'^'= 
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34000 
32000 
30000 
28000 
26000 





Ground state 




First exc. 




Second exc. 







w=2 
w=3 
w=4 - 
w=5 
w=6 



Figure 19. First three eigenvalues of the hamilto- 
nian of Harigaya's model as a function of the gen- 
eration. 



Figure 20. Optical absorption edge for extended 
dendrimers with different numbers of wedges (10 — 2 
up to It) = 6) according to Harigaya's model with 
DMRG. 



and it is obtained that Cc ~ (3190 ± 120)cm-i and oic ~ (0.54 ± 0.02). 

As a matter of fact, Harigaya's model of a single exciton is proposed only for extended den- 
drimers, although no reason is presented for which it might not be applied to compact ones. For 
the latter case, linear cluster theory predicts successfully the constancy of the data. The failure 
of the one-exciton theory leaves a single feasible possibility: the higher links concentration makes 
the effects of interaction between excitons much bigger than in the case of extended dendrimers, 
rendering the approximation incorrect. 
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4.1. Long Range Formulation of DMRG. 

In the previous chapter some "implicit" RSRG algorithms were developed for quantum me- 
chanics (QM). In these algorithms the wave- function is not explicitly stored and is not used in the 
calculations, so the RG step requires only 0(ra'^) operations, where m is the number of states we 
wish to obtain and k is a certain power. 

But the DMRG algorithm for trees -the most developed of the algorithms developed for QM- 
may not be generalized to various dimensions. This does not imply that RSRG may not be applied 
successfully to multidimensional problems, but that methods must be "explicit": the full wave- 
function must be stored and the RG step shall take, at the best possible case, 0(N) steps, with 
N the total number of sites to be considered. Even so, RSRG methods may result quite profitable 
both theoritically and numerically for these problems. 
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In Search for an Implicit Multidimensional Algorithm. 

The question which we shall succinctly review in this section is "Why is it necessary to store 
the full wave-function when the problem is 2D?" The answer is, rather briefly, because adding is 
easy, but substracting is not. The present paragraph elaborates this argument. 

DMRG is based upon the addition process or modus ponens: it was said in the former chapter 
that at each step one of the blocks grows and the other (or others) decrease, but this assertion 
should be carefully restated. Indeed, there is a block which increases size; but the information for 
the shrunken block is taken from the previous cycle. The reason is that we do not know how to 
remove sites: we do not have a substraction or modus tollens^ available. 

How does this fact influence our problem? DMRG is applicable to trees because, given any 
site, it is possible to split the system into blocks such that: 

• The blocks are disconnected among themselves. Only in this case may the modus tollens be 
simulated by taking the information from a previous step. The reason is that, in this case, the 
modification of the states of a block does not alter the matrix elements of another one, which 
conserve their validity from a sweeping cycle to the following one. 

• Each neighbour of the site lies in a different block. This must be true so as all sites may become 
"free sites" . At any RG step, all neighbouring blocks of the site but the one which contains 
the new free site shall get fused. Thus, this block should not share sites with the others. 

Once we have accepted the fact that the modus tollens or substraction process may not be 
simulated in non-tree graphs, the following question remains: "Why can not an authentic implicit 
modus tollens be built?" 

The answer to this question is the following: the removal of a site from a block forces the 
modification of matrix elements of some operators. Such a modification requires the knowledge of 
the real value of the wave-function at the site. But since the process shall be iterated, it shall be 
necessary to know, along with it, all the other values. 

The explicit algorithm is, therefore, unavoidable. It is then possible to write down a DMRG 
process for long range quantum mechanical problems, i.e.: where there is no underlying graph 
structure whatsoever. Afterwards we shall see that, if it exists, such a structure may be used to 
reduce drastically the computation time (even when it is not a tree). 

The first algorithm of this kind was developed in the work of M.A. Martm-Delgado and G. 
Sierra [MDS 99]. 

An Asymptotically Free Model in Quantum Mechanics. 

One of the physical systems which results more difficult to tackle is probably Quantum Chro- 

modynamics (QCD), which is the standard theory for the behaviour of hadronic matter (quarks 
and gluons). This theory is asymptotically free, i.e.: at short distances the interaction is so weak 
that perturbation theory may be employed. The problem appears at the energy scales at which 
bound states (hadrons) exist [DGMP 97]. The analysis is rather difficult since all scales are strongly 
coupled. 

The simplest analogue in quantum mechanics was proposed by K.G. Wilson and St.D. Glazek 



The terms modus ponens and modus tollens, inspired by the terminology of formal logic, mean in latin adding 
and removing respectively. 
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[GW 97] and it is the momentum space analysis of the behaviour of a bidimensional particle bound 
by a delta potential. 

The regularization employed for the problem required a certain IR cutoff and UV cutoff, given 
by two integer numbers M. and N. Just like the shell models of turbulence, the states have momenta 
which follow a certain geometric progression. Given a momentum scaling factor b (chosen to be 
V2), the sclf-cncrgy of the n-th state is b^"^ and the coupling between two states n and m is 
— g-v/En^m- The discretized hamiltonian therefore takes the form 

Hn.m = Snmb'^ - gb^+"^ M < Tl, m < N [1] 

The full hamiltonian has discrete scale invariance Hrt+i^Ta+i = b^Hrt^raj which is broken by the 
IR and UV cutoflFs. We say that it represents an asymptotically free model since the coupling 
between neighbour momenta are arbitrarily small in one of the extremes of the spectrum (IR or 
UV depending on the value of b). 

The analysis of Glazek and Wilson is based on the method they developed, called the Similarity 
Renormalization Group (SRG) [GW 93] [GW94I, which works by applying similarity transforma- 
tions to the hamiltonian until it becomes band diagonal. It "decouples" in momentum space the 
scales which are far one from the other in a perturbative way, by integrating Wegner's equation 



dH(s) 



[Hd(s),H(s)],H(s) 



ds 

(where H(i(s) is the diagonal part of H(s)) with the initial condition H(0) = H and then taking 
the limit s — > 00. The parameter s may be identified with the inverse square of the "width" in 
energies of the matrix H(s) and, thus, in the aforementioned limit, the hamiltonian is diagonal. 

The DMRG analysis of this problem, which is non-perturbative, makes up the first incursion 
of this technique into the field of asymptotically free theories. 



DMRG FOR Long Range Problems. 

We shall now describe the DMRG technique which was employed to deal with the given 
problem. It is a technical modification of the method introduced in [MDS 99], which fully respects 
the general idea. 

First and foremost, the system is divided into "left" and "right", which in this case means low 
momentum states and high momentum states. In any case, it is important to observe that there 
is a strong interaction between states which belong to different blocks, in contrast to what we saw 
in the previous chapter. 

Two blocks (left and right, L and R) plus two sites between them shall be considered. Block L 
is formed by the sites [1, . . . ,p — 1]. Let and {l^"*^^)} with i G [1, . . . ,m] be the series of 

m orthonormal states for the left and right sides. At each step the full wave-functions should be 
available: 

{(5p|tl;">} {(5q|xl;'^^>} 
for all p e [1, . . . ,p — 1] and q € [p + 2, . . . , N]. I.e.: in total m(N — 2) numbers. 

The Ansatz is exactly the same as in the ID-DMRG case: 



86 



4. Multidimensional Formulations of the DMRG. 



let)) = x; i^"> + c™+' i5p) + i6p+i) + x; i^''') 



1=1 



1=1 



With this Ansatz we may build a superblock hamiltonian which is different from the one we 
previously used: 



H 



Sb 



v 



L 




yLrt 


HLR 


jLl 


Hcl,cl 








^CT,cl 


^^CT.CT 


jRr 




jRlt 


jRrt 


R 



[2] 



In this hamiltonian the classical elements of the DMRG are present, which are: 



• The intra-block elements = H |^^''-'), Rij = H 

• The hooks of the blocks to the free sites: T^^ = H |6p), T^^ 

• The terms corresponding to the free sites: Hcv,cl) etc. 



But there appear three new elements of the superblock matrix which are non-zero a priori: 

^il)*^^! H |5p+i) and viceversa 



Lr 



• The link between the left part and the right free site 
TRi = (xl;Ri|H|6p). 

• The link between the left and right parts: H[--'^ = (i|>'-^| H |it>'^'). 

All these matrix elements of the superblock hamiltonian may be fully calculated if the states 
are available, but it is convenient to hold in memory the superblock hamiltonian of the former step 
and modify it consequently. 

The superblock hamiltonian is diagonalized and its m lowest energy states are retained. We 
now give two options for the last step of the procedure, which is the states updating. The first is 
more simple but much less efficient. In both cases, we shall assume that the left block is growing. 

• First procedure Let us build the complete states \^^) by substituting the coefficients of the 
eigcnstates of H^'' into the Ansatz. We write down these m states twice. On the first copy we 
make all components from the p + 1 onwards vanish, while on the second one we put zeroes on all 
components up to the number p + 2. 

These states correspond to the new L and R states respectively, but they require re-orthonor- 
malization through a Gram-Schmidt process. Once it is finished, the matrix elements of the new 
superblock hamiltonian H^'' are computed. 

The full process requires O(N^) operations. 

• Second procedure. A better notation is required to explain this method, which is much more 
complex (and efficient). Let capital indices 1, J, . . . run over the left half of the components of the 
eigenstates of H^'', in the range [1, . . . ,m + 1]. Let us define the original states of the left part 
(with "Block + •" structure) as 
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ULiN n^l^^^) ifl<m 
1^ / \ |6p) ifl = m+l 

I.e.: we include the state |6p) into this set. Let us consider now the basis change matrix obtained 
through the eigenstates of H^'^: Bj denotes the J-th component of the i-th eigenstate (with i G 
[1, . . . , m] and J e [1, . . . , Ta+ 1]) . By not allowing J to run over its full natural range (which would 
be [1, . . . , 2m + 2]), we are losing the orthogonality of the transformation. Let us see how we may 
recover it. 

Let Cij be the dot products matrix among the states of B, defined by 

m+l 
K=l 

Appendix D shows how to obtain the Gram-Schmidt basis change rather fast by using this 
information. Let G be the matrix implementing that basis change and let us define 

Bp = GjB} 

to be the orthogonal basis change matrix from the old states \ ^^^) to the new ones, given by 

m 

3 = 1 

It is possible to write an equivalent basis change for the right side states (with structure 
"Block — •", i.e.: block minus site). But there is an added difficulty: the state to remove (|6p+2)) 
is not one of the old states of the basis of H^''. Let us see how to approach the new problem. 

Let 4'p+2 be the component of the i-th state of the right block on site p -|- 2. Such a block 
undergoes the transformation 

which we shall call substraction algorithm (or modus tollens). These states arc not orthonormal. 
So as they can be, the fast Gram-Schmidt procedure is applied (see appendix D) on its dot products 
matrix: 

and let us call Q| the basis change matrix which ensures orthonormality. Defining the right side 
states as 



\^\>^^) ifl<m 
|6p+2) if I = m -1-1 



and the transformation matrix 



Bi< = 6i<-t^pV26K.m+l 

the new right side states may be expressed as 
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Therefore, we have formally the same type of transformation for both blocks, even though we 
know that there are rather different operations involved. The following step is the computation of 
the new matrix elements between those renormalized states. Let us see how the calculation may 
be done for each box the matrix H^'^. 

• Intra-block left part. The computation is straightforward. 

m+l 



Li 



H 



1,1=1 



But since all the states appearing in the sum are states of the old H^*^, we have only to transform 
the upper-left corner of the old matrix: 

m+l 
ij=i 

• Left hook to the left free site. The new free left site is p + 1, which beforehand was the right 
free site. Thus, 



T 



Ll 



Li 



m+l 



H|6p+i)= 5^Bp<cl)"|H|6p+i) 



1=1 



In the same fashion as before we obtain 



m+l 



m+2 



1=1 



• Link between the left block and the right free site. The new right free site did not exist before 
as an independent entity, so there are no shortcuts for the computation: 



-Lr 



^1^ 



Li 



H|6p+2) 



The only possible "trick" to reduce the computational effort comes up if the hamiltonian 
respects some kind of neighborhood structure (even if it is not linear). 

• Left-right inter-blocks part. The calculation provides some terms which may be simplified 
and other ones which may not. 



H 



LR 



Li 



H 



m+l 



1,1=1 



E Bt^Bj:' {^^^\ H + E ^^+1 (^'1 ^ |Sp+2> + E BJnViBr (6p| H |xl;^^> + 

k,l=l lc=l 1=1 



+Bj,ViBS+i(6p|H|6p+2) 
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Fortunately, only the second term requires a full computation (because it makes interact the 
site p + 2 with the left block). Thus, we have 

mm m 

lc,l=l k=l k=l 

• Right intra-blocks part. This part has also got some term which must be calculated in full 
(of course, it involves site p + 2). 



Ri 



H 



m 



k,l=l 



m 

+ 5^ Bf B^+i (il;'^'^! H |6p+2> + (i ^ i, he) + B^ViB^+i (6p+2| H |6p+2) 
k=l 

m m 

Y: B^Bf'Ru + E (b^bS+i + B-^^B^V,) Re {V'^l H \5^^,) + B^V,BS+iHp+2.p+2 
c,l=l k=l 



• Right link to the left free site. This is straightforward, since the site |6p+i) belonged to the 
old superblock. 



Rl 



^1^ 



Ri 



H |6p+i) = £ <-i|;'^''| H |6p+i) + B^Vi (5p+2| H |6p+i) 



k=l 



°\c 'k + t>m+l'^P+2,p+l 



k=l 



• Right hook to the right free site. This new free site is, of course, |6p+2), so we must compute 
fully some matrix elements: 



H |6p+2) = X] {^^^\ ^ |Sp+2) + Bj^\iHp+2,p+2 



k=l 



We may observe that the computations we have named as full only require the knowledge of 
2m numbers: 



(xl;'-^|H|5p+2) (i^''^|H|&p+2) 
Thus, the calculation of the matrix elements (and the updating of the wave-functions) is 
performed in 0(N) steps, while the naive version was carried out in 0(N'^) operations'^. It is 
the computation of the second matrix element \ H |6p+2) which forces the storage of the full 
wave-functions. 



Notice that in the short range ID-DMRG the first of these matrix elements is null, while the second one 
would have been taken from a previous cycle, leaving the RG step with 0{m?) operations. 
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Application to the Asymptotically Free Model. 

The implementation of the original long range DMRG to the 2D delta potential in momentum 
space (considered to be a quantum-mechanical model of an asymptotically free system), which 
not only lacks an underlying graph structure, but also combines matrix elements of very different 
orders of magnitude, had a rather good precision. 

As a test, system [1] was considered with n = 38 "shells"^ with a scale factor b = \/2, IR 
cutoff M. = —21 and UV cutoff N = 16. The coupling constant is fixed to a precise value^ so as to 
obtain a ground state with energy —1. The obtained wave-function appears in figure 1. 




Figure 1. Ground state of the asymptotically free model under study. The abscissa represents the momentum 
"shell". 

As a result it yields, in two sweeping cycles, a precision of 15 orders of magnitude for the 
ground state and, if m = 4, 6 orders of magnitude for the three first excited states. For more 
details, calculated with the first version of the long range DMRG (which gives exactly the same 
results albeit more slowly), see [MDS 99]. 

Application to Short Range 2D Hamiltonians. 

The case of a bidimcnsional problem with short range terms is not analyzable by the DMRG 
which was developed in the former chapter, since the underlying graph has no tree structure. 
With the formalism exposed in this section, it is possible to study such systems with an arbitrary 
potential, and we may take profit from the underlying graph structure so as to compute quickly 
the necessary matrix elements of the hamiltonian. 

The method starts with a "unidimensionalization" of the bidimensional system, converted into 
a snake, as it is shown in figure 2. 

The links which, in a bidimensional lattice, are local and natural, appear under this transfor- 
mation to be long ranged. 2D-DMRG computations for many-body problems yield a much lower 
precision than in the ID case [WHI99]. The left-right distinction, which forces unidimensionaliza- 
tion, is unnatural and destroys the good properties of DMRG. 

The reader may ask herself whether n = 38 is a respectable number. The answer is that, in these models 
where the shell momenta follow a geometric progression, one runs the risk of overloading the machine precision. 
38 jumps of a factor b = \f2 is already a very interesting computation. 

^ 5 = 0.0606060003210886 
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Figure 2. Unidimensionalization of a 3 x 3 bidimensional lattice. The "snake" gives the order of the sites in 
the new structure, and the broken links appear in dotted Hues. On the right, the hnear chain with its long range 
structure. 

The numerical results for this technique shall appear in section 4.4. The convergence of these 
results is (at least) as safe as those in ID, but the computation times are bigger. 



4.2. Punctures Renormalization Group. 

The long range DMRG algorithms are based on the left-right distinction, which only makes 
sense in ID. This distinction is unnecessary and even harmful in many cases, as it shall be clear 
when comparing the results of the 2D DMRG with the new method which is about to be introduced 
in this section. 

This new method, which has been called Punctures Renormalization Group (PRG), was de- 
veloped by M.A. Martm-Delgado, G. Sierra and the present author [MRSoo]. It requires a single 
block instead of two. Some "punctures" are drilled into this unique block, which allow us to study 
with detail the behaviour of the system at individual sites. 

In the first section the simplest version of the algorithm is explained, in which a single state 
and a single puncture are studied. 

Single Block - Single State - Single Puncture. 

Let us consider a global state \^\>o), defined over the sites of a given graph G- This state is a 
certain approximation to the ground state of a hamiltonian which fulfills the connectivity rules of 
Q. We may obtain a better approximation to the ground state by following this process: 

• We choose any site p of the graph. 

• We project on the subspace orthogonal to |5p). Or, in other words, we make the component 
\|jo(p) vanish. We re-normalize the state. Let us call it now |i|Jo)- 

• We write an Ansatz of the following form: 



• As the states entering the summation are orthogonal, we may proceed to write an effective 
hamiltonian (for the superblock): 





) 
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• We diagonalize the hamiltonian and retain only the ground state. This has only two compo- 
nents: the weight of the block and the weight of the state |6p). The lowest eigenvalue of that 
matrix is an approximation to the ground state energy. 



It is convenient to emply a graphical representation. Let us consider a bidimensional lattice 
such as that of figure 3. 



Figure 3. Representation of the "single block + single puncture" scheme. 

Thus, the block contains the whole system but a site, which we shall call the puncture. Figure 
3 shows that the block (shadowed) is connected, so it makes no sense to split it. This is the "single 
block" philosophy, essential to the method. 

When we diagonalize the superblock hamiltonian we are looking for the best approximation 
of the ground state (and, incidentally, of the first excited state) in the subspace spanned by the 
punctured block and the delta state at the site. An appropriate way of thinking about the process 
is that the wavefunction improves at the site, meanwhile the rest of the system is merely adapted 
to the norm changes which are introduced. 

After a RG step is finished, the puncture may translate to any other site of the system. The 
best results, notwithstanding, are obtained when the sites through which the "puncture" travels 
make up a connected path on the graph. Moreover, let us state that the movement from a site p 
to a neighbour q covers the link (pq); it is then convenient to cover all the links of the graph in 
at least one sense. Let us call a PRG- "sweeping" to a sequence of steps such as the former one 
which passes through all the sites, leaving the covering of all the links as a desideratum. If this 
last requirement is fulfilled we shall call it a "sewing" . 

Single Block - Various States - Various Punctures. 

Once the basic idea of the method has been exposed, we give the technical details for the 

process with many states and punctures, such as it is performed in practice. The system takes 
profit of any graph structure, but it is also efficient under its absence (as in the case of the delta 
potential formerly described). 



Any orthonormal set of states is suitable for the beginning of the computation. Obviously, 
a right guess speeds the convergence up, as it is always true for the variational methods. We 
shall discuss in section 4.5 various appropriate warmup techniques, but for the moment we shall 
undertake the most simple one: the orthonormalization of a set of states whose components are 
randomly selected. 
• Punctures and Sewing. 

The set of punctures (or patch) will be denoted by P. Its size, Np, may change throughout 
the process. The strategy for the election of a patch and its movement through the system is a 




• o o o o 



• o o o e 




• Warmup. 
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subject which, inasmuch as the warmup, may only affect the convergence speed. 

Just so as not to leave everything living in an abstract realm, let us propose a concrete example 
for ID. At the t-th PRG step, the patch may be given by the set P(t) = {2t,2t + 1}. Thus, the 
patches would all be made up of two punctures. Notice that the set P(t) = {t, t + 1} is very similar 
but not equal to the former one. In the second case, patches overlap. Therefore, more PRG steps 
are required in order to complete a sewing. It may be worth if the number of necessary sewings is 
substantially smaller. 

For 2D lattices a customary patch is the 2x2 square: 

P = {(^J).(i+iJ).(i.i + i).(^ + i.i + i)} 

All the example patches which have been given are formed by connected punctures, although 
it is not a sine qua non condition. 

• States, Ansatz and Superblock. 

Once the patch P has been chosen, the Np delta states {|5p)}pGP become part of the Ansatz. 
The block states shall be obtained from the states which were achieved in the last step, denoted 
by {l^^i)}!!^!- Of course, m is the number of states we wish to obtain. 

Let us call Qp the projector on the states which are orthogonal to the set of delta states on 
the patch sites: 

Qp = I-5^|5p)(5p| 

The action of this operator is easy to describe. It simply removes the components which 
correspond to the patch sites: (6p| Qp = for all the p G P. 

Maybe the most important technical problem of the procedure appears at this point: the 
states {Qp lij^i)} s-re not an orthogonal set. The re-orthonormalization may be carried out by 
the standard Gram-Schmidt technique, but it is more appropriate to use the fast Gram-Schmidt 
procedure described in appendix D. It is thus important to know the scalar products matrix 
Cij = (i^il Qp li^j). Let G-- be the basis-change matrix obtained through the process. Now the 
states 

m. 

J]g|Qp|i|;0 ie[l,...,m] 

j=i 

make up an orthonormal set. These states, along with the delta states of the patch, make up the 
Ansatz: 

m. Np 

14)) = 5^ ai 111;?) + a^+^ |6p(j)) [3] 
t=i j=i 

where P(j) denote each of the patch punctures. The {cii}^^'^^ are the variational parameters. 

The superblock hamiltonian Hsb which corresponds to the Ansatz [3] is now written straight- 
forwardly: 
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/ /i^liili^V 

/ (i|)l|H|il;^) • 


• (^l>m|H|tl)t) 


1 c 1 1 ■ 1 1 ^ V 

(6p(i)|H|tl)t) • 


• <5p(Np)|H|i|;t) \ 


(M^JIHiti;;,) • 


• (^m|H|^l);,) 


<6p(i)|H|ib;,) • 


• (5p(Np)|H|il;^;,) 


(il)*|H|6p(i)) • 


• (^^;^|H|6p(l)) 


(Sp(i)l H |6p(i)) • 


■ (^P(Np)l ^ I^P(l)) 


V (tl)t|H|6p(N^)> . 


• (MJm|H|6p(N^)) 


(5p(i)| H |6p(Np)) • 





[4] 



The structure of the matrix [4] is quite clear: 



H 



Sb 



H 



t 

BP 



H 



BP 



where Hb is the single block hamiltonian (the |iJj*) among themselves). Hp is the set of elements 
related to the patch, directly taken from the total hamiltonian. Hbp are the mixed elements, 
between block and patch states. 

The aforementioned matrix elements may, of course, be computed ah initio^ i.e.: from the real 
states components. But, again, this involves an inneccssary waste of resources since: 

• We know the superblock hamiltonian at the former PRG step, and it may be "merely" adapted 
to the new circumstances. 

• An underlying graph structure boosts the computation of the matrix elements. 

If it is naively computed, Hsb requires (m + Np)'^N'^ operations. We shall assume, therefore, 
that the former step Hsb is known: 

Hi,- = (ibi|H|il;j) 

From this expression^ we obtain Hsb in the straightest possible way. Foremost, we define the 
intermediate matrix 

^ m QpHQp \M) 

Of course, these numbers are not matrix elements of the hamiltonian for any set of real states, 
since the {Qp |i]^t)} are not orthonormalized. Nonetheless, they are easy to find: 



= ^i(i<)Hicix^i(i) = E-E-E+E 

k,leS-P I k.leS keP kes 1c,leP 

\ les lep 



where S denotes again the set of sites of the system. The first summation runs over the sites of the 
block. Using the simplest version of the ^^inclusion-exclusion theorem" , it is split into four easier 
partial sums: 



For a general PRG step, hij = EiSij, i.e.: matrix h is diagonal and its elements are the estimates for the 
energy at the former step. Nevertheless, matrix h may be "complete" after the warmup, so we conserve the 
general form. 
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• The first summation is, of course, Hij . 

• The fourth one is rather fast to compute, since it only involves matrix elements among states 
of P. The number of operations is 0(Np). 

• The second and third summations, equal in structure, must be computed in full if there is no 
underlying graph structure. Otherwise, they are fast to calculate: 

5^il;i(k)HuX^j(l) = x|;t(k)Huil^j(l) 

where N(k) denotes the set of neighbours of k. These summations require 0(z x Np) operations, 
where z is an upper bound for the size of N(p) for all p (maximum coordination index). 

After the computation of h,*^ , it is easy to obtain Hb : 

m m 

k,l=l k,l,= l 

The elements of Hbp are also quick to compute under a neighbourhood structure: 



m m 

(HBp)ij = m H |6p(j)> = J2 Gic (^kl QpH |5p(3)> = Y E Gi,tl)„(s)H3p(„ 

k=l k=l s€N(P(j))-P 

where the last summation over N(P(j)) — P means "the neighbours of the j-th puncture which are 
not punctures themselves" . This computation takes also less than z x m operations. 

After all these operations, all the elements of Hsb have been computed without using in any 
case all the components of the states. The neighbourhood structure reduces the timing drastically. 
If it is absent, the computation takes 0(N) operations (again, N S> m is assumed). 

• Truncation and States Recomposition. 

The superblock hamiltonian represents the whole of the system. When it is diagonalized, its m 
lowest energy states provide us with the variational parameters {ci\} (with j £ [1, . . . , m] -denoting 
which excited state- and i G [1, . . . , m + Np] -denoting which component) which must be inserted 
into [3]. These values yield the best approximation, within the subspace spanned by our states, to 
the lowest energy states of the system. The eigenvalues of Hsb give, as it was assumed, the best 
estimates for the lowest energies of the full system. 

Using these values we reconstruct the wave-functions: 

ra Np 
i=l }=1 

But the stored states are the |Tjji) and not the \^\>l)- Therefore, 

4'k)= E 4G>Qp|xl;j) + Ea^+j|5pa)> [5] 
ij=i j=i 
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The action of the operator Qp is really simple: make all the components associated with P 
vanish. We may build just a matrix which encapsulates the re-orthogonalization and the renor- 
malization: 



i=l 



where the matrix G* is defined by 



[ G- if i, j < m 
(G*)^ = < 6} if ij >m 
. otherwise 
The reconstruction procedure yields, therefore: 



And, under a neighbourhood structure, this is the only step which scales with N. Thus, in the 
worst case scenario, the number of operations for the PRG step is 0(N). 

There is only one more detail to be closed: the matrix 

H 



this 



is necessary to consider at the next step. Due to the construction rule of the states 
matrix is diagonal, and its elements are the energies of the states: 

y^.^fU ifi = i 

1 otherwise 

We have mantained the notation Hij for the sake of generality (e.g.: the states coming out of the 
warmup need not have the same structure). 

Summary. 

• Stored Data. The N components of the m states li^ii) and the elements of hij (which usually 
make up a diagonal matrix). 

• System States and Ansatz. The punctures components are removed with the operator Qp. 
The states are re-orthonormalized with a Gram-Schmidt matrix G. New states are defined 



l*i) = 

The final states are created after the Ansatz 



1=1 1 

6p(i_ra)) ifi>m 



m+Np 

14))= 

i=l 

where the {ai} with i G [1, . . . , m + Np] are the variational parameters. 
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• Superblock Hamiltonian Construction. Given by Hsb = (cf'il H|4)j). Its construction, using 
the neighbourhood structure and the knowledge of h|j is rather fast. 

• Diagonalization of the Superblock and Updating. The m lowest eigenstates of the superblock 
hamiltonian are retained and the at weights are inserted into the Ansatz according to the expression 
[5]. The diagonal matrix elements h,|j are the obtained energies. 

General Features of the PRG. 

The Real Space Renormalization Group method known as PRG {Punctures Renormalization 
Group) has the following general features: 

• Ensured Convergence. Section 4.5 discusses the conditions under which such statement may 
be held. 

• Absolute Generality. Any quantum-mechanical model in a finite dimensional Hilbert space is 
analyzable through this technique. 

• Adaptative Efficience. It takes profit if the hamiltonian matrix is sparse, i.e.: the existence of 
an underlying neighbourhood structure. 

• Explicit Method. Wave-functions are stored in full. The number of operations of a PRG step 
may be not lower than 0(N), the inequality being saturated in the case of underlying graph 
structure. 



4.3. ID Implementation and Numerical Results. 

For the unidimensional examples we have chosen the following models: 

• Free particle in a box divided into N = 200 cells with fixed boundary conditions. The eigen- 
states are extended throughout the system and the boundary conditions are rather relevant. 

• Particle in a harmonic potential, with space again divided into N = 200 cells. The boundary 
conditions are much less relevant and the eigenstates are spatially localized. 

• Free particle in a ID box with periodic boundary conditions. The underlying graph is not a 
tree. 

In all the cases the objective of the calculations were the ra = 4 lowest energy eigenstates and 
the convergence criterion was a precision of one part in 10-^° for the two first states compared to 
the exact values^. 

Table 1 shows the numerical results. The numbers given in parenthesis after the PRG heading 
are the number of punctures for the patch. Respectively, 4 and 10 punctures were used. 

The numerical results for ID yield a rather clear advantage of DMRG with respect to the other 
methods on graphs "without cycles". The reason is obvious: DMRG is implicit and, in ID, the 
partition into left and right blocks is natural. PRG in ID gathers in a single block two regions 
which are really separated in space, slowing down the convergence. 

It must be remarked the case of periodic boundary conditions (third column). The underlying 
graph structure is not a tree for this system, since it contains a cycle. Therefore, the left-right 



Obviously, there are also convergence criteria which are fully internal, as it is the imposition that the variation 
through a full sweep to be smaller than a prefixed quantity. 
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Fixed B.C. 



Harm. Pot. 



Periodic B.C. 



Exact 
DMRG 

PRC (4) 

PRC (10) 



0.83 (3) 
25.6 (198) 
6.97 (65) 



2.13 



0.6 (2) 
9.8 (76) 
3.22 (30) 



2.14 



109 (11) 
10.94 (87) 
3.81 (34) 



2.09 



Table 1. Benchmarking for the ID-PRG algorithm, DMRG and exact diagonalization for a) a free particle with 
fixed b.c, N = 200 and m = 4, 6^ in presence of an harmonic potential, c) free, but with periodic boundary 
conditions. In all the cases a precision was reached of one part in lO"^*^ for the first two states (when the exact 
energy was null, the precision was imposed in absolute value). The numbers show the CPU time in arbitrary 
units (Pentium III at 450 MHz) and, in parenthesis, the number of PRG or DMRG sweeps. 

distinction is unnatural and DMRG leads to worse results than PRG. 

Both methods take a worse result than the exact diagonalization for N = 200, but for bigger 
sizes (empirically up to N = 120.000 sites), due to the different scaling regimes, PRG becomes a 
rather suitable method. CPU times scale for all methods (exact, DMRG and PRG) as a power 
law: 



Albeit the parameter K may be rather relevant for practical applications, it is usual to consider 
the exponent oc to be the key, due to its "universality" . In other terms: an improvement in the 
implementation or in hardware may lower K, but the exponent a may only diminish with a deep 
change in the computation algorithm. This statement is illustrated by the following results: 

• Exact diagonalization (periodic b.c): a 3.1. 

• Implicit DMRG (fixed b.c): oi^ 1.2. 

• Explicit DMRG (periodic b.c): a 3.7. 

• PRG (periodic b.c): a ^ 2.2. 

It may be observed that the explicit DMRG for non-tree graphs scales worse than the exact 
diagonalization. DMRG in its own scope is, no doubt, the one which presents the best performance. 
PRG is an alternative when DMRG can not tackle the problem in an appropriate way: for long 
range and/or multidimensional systems. 

Some practical questions on implementation: 

• When the number of punctures does not divide exactly the lattice size, we may choose between 
these options for the analysis of the last patch, a) The last patch contains less punctures than 
the rest: Np = mod Np, ot b) the movement backwards starts when there are less than Np 
sites left to reach the extreme. 

• The movement of the patch may be carried out with or without punctures overlapping. In 
ID the first case appears to be slightly more efficient for big patches (r^ 10). Figure 4 shows 
graphically the meaning of advance with overlapping. 



tcpu ~ K • N 



a 



4.3. ID Implementation and Numerical Results. 
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Figure 4. Movement of a 4 sites patch with overlapping in a ID lattice. On the first PRG step, the patch 
includes sites 1 upto 4. On the second step, sites 4 to 8. There is a common site between successive patches. 



Convergence in DMRG and PRG. 



Focusing on the DMRG and PRG techniques for ID systems {without periodic b.c), figures 
5 and 6 show the different approaches to convergence for two free systems with fixed b.c. and 
N = 100 sites. 
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Figure 5. DMRG Convergence. Relative error for the 
first four states as a function of the RG step. Notice 
how convergence is led by fast downwards slopes and 
long plateaux. 



Figure 6. PRG Convergence. The logarithm of the 
error is presented as a function of the RG step. The 
straight lines show the fit to an exponential decay. 



DMRG advances by very fast strides of the error and rather long plateaux. On the other hand, 
PRG is a "long-distance runner" , as it is shown by the good fit to a law of the type 

AE » exp(-K • n) 

(where n is the RG step index) in the whole range of values until machine precision is reached. 

This difference may be explained in qualitative terms. DMRG uses for the decreasing block 
the information of a former step. This implies that at the start of a "half sweep" (e.g. a sweep 
left — > right) we are employing a very big block which has not been recently updated. The biggest 
advances are made when the decreasing block is very small (in figure 5, the high slopes are just 
before steps 200 and 400, which correspond to a sense change). 

On the other hand, PRG has all its information updated. At every step, the change in the wave- 
functions is proportional to its distance to the desired state. This explains, at least qualitatively, 
the exponential convergence. 

It is interesting to study the PRG "residual" state. Let us consider the case when only the 
ground state is retained (m =1), whose exact value we shall assume to be known. We substract 
at each RG step the exact state from the approximate one and normalize the result. This residual 
state, after a number of RG steps, becomes the first excited state. The transient is longer when 
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both states have different symmetries (even vs. odd). In this case, the residual state corresponds 
to the second excited state for some time. 

When m states are conserved, the residual of the i-th state converges to the (ra + l)-th state 
after a transient, which may be retarded for symmetry reasons (even vs. odd). 



4.4. PRG Analysis for 2D and 3D Lattices. 

The title for this chapter reminds us that the fundamental reason for the improvements on 
the DMRG was the study of 2D and 3D systems. This section deals with a series of fundamental 
technical questions when applying the PRG in 2D and 3D, along with a series of applications and 
numerical results. 



Sewing Lattices. 

Let us consider a 2D rectangular lattice x L^ and a puncture formed by a single site. For 
DMRG the only possible path is the one we shall call "DMRG sweep" , which traverses the system 
in an anisotropic zig-zag (see figure 2). 

PRG, dealing with a single block, allows an isotropic path in which not only all sites are visited, 
but also all links. In the terminology introduced in section 4.2 we may talk of a "sewing". 

The change from the anisotropic to the isotropic sweeping obeys the following reason: the 
derivatives of the wave--functions along a link between two sites are updated more appropriately 
when such a link is traversed by the puncture. The DMRG sweeping leaves almost half of the links 
without covering, slowing down the convergence. 

The algorithm for drawing a sewing depends on the parity of L^ and Ly. Figure 7 shows a 
sewing for the even-even lattice. 



Figure 7. A possible sewing for a 4 x 4 lattice. 



The structure of the sewings is exposed more rigorously in an algebraic notation, which on the 
other hand describes the code of the employed programs. Let us define the four operators L, R, U 
and D as the moves of the puncture leftwards, rightwards, upwards and downwards respectively. 
The path in figure 7 is algebraically described in pseudocode 1. 

The odd-odd sewing, exemplified in figure 8 for a 5 x 5 lattice, is described in pseudocode 2 in 
a more rigorous way^. 

^ The last expression of pseudocode 1, describing the sewing of an even-even lattice in terms of two operators 
and their inverses, reminds the topological description of a torus. On the other hand, the description of the 
odd-odd sewing (pseudocode 2) would be topologically trivial [BLA82]. 



4.4. PRG Analysis for 2D and 3D Lattices. 
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Pseudocode 1. Algorithm for the sewing of an even-even lattice. 



Figure 8. Sewing for a 5 x 5 lattice. Notice that, in this case, each link is traversed in both senses. 



[,^uLy-lRuLy-l](L,-l)/2 

L'-"-! [UR^^-^UV-^-^]^^y-^^/^ 
7-1 7-1 7 7 

■^hor ■^vert ■^vert ■^hor 

Pseudocode 2. Algorithm for the sewing of an odd-odd lattice. 

As a matter of fact, the sewing is one of the most important features of the bidimensional and 
tridimensional PRG. 2D DMRG, due to its own unidirectional structure, is forced to carry out the 
sweeping in a necessarily anisotropic way. 

Results improve significantly when patches composed of various punctures are employed. Hav- 
ing tried many alternatives, we have checked that the square patch of Lp x Lp sites is the one 
leading to best results. The algorithm exposed previously for the sewing is suitable when Lp is 
commensurate both with Lx and Ly . In this case, it is possible to split the system into blocks of 
size Lp X Lp for the patch to cover. Figure 9 illustrates the process. 




Path 



Figure 9. A path through a 6 x 6 lattice by a 2 x 2 patch. The dimensions of the patch must be commensurate 
with the ones of the lattice. 
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Numerical Results in 2D. 

Two quantum-mechanical problems have been chosen so as to check the performance of the 
2D-PRG: the free particle in a box (with fixed b.c.) and the bidimensional hydrogen atom. 

The first one is fully analogous to the ID problem: the diagonalization of the laplacian matrix 
of the 2D graph with fixed b.c. The eigenstates are analytically obtainable. If = Ly, the D4 
group of symmetries guarantees the degeneration of the first excited state (among many others): 
there are two twin states with nodes lines which may point in any couple of orthogonal directions. 

The hydrogen atom on a 2D lattice is exactly solvable in the continuum limit. If and 
Ly are not both odd, then the atomic nucleus stays at the center of a "plaquette" or an "edge" . 
Otherwise, a regularization is required so as the potential energy is finite throughout the system. 
In the continuum limit this distinction is senseless and energies are 

E^^ = -ZVn^ n = l,2,... 

The production of a graphical output in the form of a "movie" is straightforward within PRG. 
Figure 10 shows a typical picture from that film. 




Figure 10. A density plot for the wave-functions obtained with the 2D-PRG algorithm for the particle on a 
15 X 15 lattice under the effect of a Coulomb potential, before convergence has been reached. The states grow in 
energy in the sense of reading. White pixels at the left side of the eight states correspond to the position of the 
2x2 patch. 

The computation time has been obtained for the hydrogen atom through three methods: 2D- 
DMRG, PRG and exact diagonalization. The results for different lattices are shown in figure 11. 
The three data series fit to power laws: 

tcpu fx L*^ 

Figure 11 shows that the lowest slope corresponds to PRG. In table 2 we may observe the 
corresponding scaling exponents. These results prove that PRG is the best choice for big lattices, 
since its convergence is assured and it is faster than the others. 



4.4. PRG Analysis for 2D and 3D Lattices. 
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Figure 11. Log-log graph of the CPU time versus the length L of a square lattice in a Coulomb potential. Four 
states are stored, a precision of 10"^ is required for 2D-DMRG and 10"^° for PRG with a 2 x 2 patch. 





Free Part. 


2D Hydrogen 


Exact 


6.4 ±0.3 


6.61 ±0.27 


DMRG 


5.5 ±0.1 


5.21 ±0.05 


PRG 


3.92 ±0.1 


3.85 ±0.07 



Table 2. Scaling exponents for the free particle with fixed b.c. and for the particle in a Coulomb potential. 

PRG IN 3 Dimensions. 

The PRG algorithm is generalizable to any dimension. Even more, results improve compar- 
atively with it. We have prepared programs for a 3D lattice, although tests have not been so 
exhaustive. 

The 2D sewing algorithms have their 3D analogues, which are naturally more complex. Thus, 
for the sake of simplicity, we shall only show the path for the odd-odd-odd case^, shown in figure 
12 and pseudocode 3. 

Pseudocode 3 uses the following movement operators: B and F (x-axis), L and R (y-axis), D 
and U (z-axis). Values of the exponent 9 for the hydrogen atom are given in table 3. 



Method 


9 


Exact 


9.5 ±0.6 


PRG (2) 


6.6 ± 0.4 


PRG (3) 


5.6 ± 0.3 



Table 3. Scaling exponent 6 for the CPU time versus the linear size of the 3D lattice (tcpu oc L ). The number 
in parenthesis shows the cubic patch size. 



Such a solution was initially developed by Silvia N. Santalla, who put it kindly at our disposal. 
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Figure 12. Three-dimensional sewing for an odd-odd-odd lattice, exemplified with a 3 x 3 x 3 one. Red lines 
are parallel to the a;-axis, green ones to the j/-axis and blue ones to the 2-axis. Notice that each link is traversed 
at least once, and inner links are covered just once. The notation Px^ and Qx^iPx ^) is defined in pseudocode 3. 



4.5. Warmups for the PRC. 



Kadanoff Partition into Blocks. 



Although it does not affect the final convergence of the RSRG explicit algorithms (long range 
DMRG and PRC), a good warmup may speed up the process considerably. A technique which may 
provide the basis of a full RSRG algorithm on its own (see following chapter) is the partition into 
blocks a la Kadanoff, also known as "Kadanoff blocking" or "Kadanoff coarse-graining" , which we 
shall describe in this section. 

The simplest presentation of this technique is for the ID system with N sites and fixed or free 
b.c. Let H be the hamiltonian N x N matrix and M be any divisor of N. Then, the complete 
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define 


Z(P,Q,a,b) = pa-l[Q(p-l)a-lQpa-l](b-l)/2 




define 


Qz = (B^L;F^R) 




define 


= (R^U;L^D) 




define 


Qy = (D^B;U^F) 




Pz ^ 


Z(R,F,Ly,U) [U Z(B,L,Lx,^) U Z(R, F, Ly , 


l)/2 


Px ^ 


Z(D,L,U,Ly) [B Z(R,U,Ly,U) B Z(D, L, U, ^)](Lx 


-l)/2 


Pv - 


Z(F,U,Lx,U) [R Z(D,B,U,Lx) R Z(F,U,Lx,U)](^^- 


-l)/2 


Path ^ 


- Qv(Py Qx(Px') Qz(Pr') Py Px Pz 





Pseudocode 3. Algorithm for the sewing of an odd-odd-odd lattice shown in figure 12. Abbreviations L — R, 
D — U and F — B correspond to movements in each axis {left-right, down-up, forward-backward). Operators Qx, 
Qy and Qz denote appropriate rotations around each axis which swaps the direction operators. Obviously, it is 
accepted that = D, = L and = F. Between the path Px^ and QxiiPxj^) a sewing of all planes 
perpendicular to axis Xi is carried out. 



Figure 13. Making blocks within a iV = 20 sites ID graph. 

system may be split into M cells of equal size f = N/M. 

A characteristic wave-function |xt) for each block is defined, with I G [1, . . . M]. 

Ixi> = ^7 E \h) 

In other words: normaHzed states which are uniform over each of the blocks. After that, the 
effective hamiltonian is 

Hfj = (Xi|H|xi) 

which is a M X M matrix. The computation of the matrix elements is immediate. Let us define 
the sets Si = [(i — l)f + 1, . . . ,lf], for which every pair is disjoint and whose union makes up the 
full lattice. Thus, 



H5 = iE E"' 



If the original hamiltonian H is the laplacian on the graph, the effective hamiltonian H*^ is 
called the collapsed laplacian in graph theory [BOL 98]. 

The eigenstatcs of H*^ arc posteriorly expanded to yield states for the full graph. Let a-' be 
the j-th component (j G [1, . . . M]) of the k-th eigenstate of H*^. Then, the full states are 
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It is interesting to ask about the similarity between the functions so obtained and the exact 
eigenstates of the hamiltonian. Let us see an example in figure 14: the free particle in a box 
discretized into N = 108 sites, divided into M = 6 blocks. The number of states to be considered 
is m = 4. 





Figure 14. Results of the partition of a ID lattice with 108 sites split into 6 Kadanoff blocks. The smooth 
curves represent the eigenstates. 



It may be observed in figure 14 that the "coarse-grained" wave-functions resemble the au- 
thentic ones as much as it is possible within their own subspace, which is a consequence of the 
variational nature of the computation. Nevertheless, the energies of the Kadanoff states and the 
exact ones differ enormously, as it may be checked in table 4. 





Ground state 


1^* exc. 


2^^! exc. 


3""^ exc. 


Kadanoff 


0.0110035 


0.0418345 


0.0863866 


0.135836 


Exact 


0.000830647 


0.0033219 


0.00747169 


0.0132766 



Table 4. Energies of the wave-functions in figure 14. 

The similarity between the warmup and the exact states may be measured through the £? 
norm of its difference. We obtain for the states of figure 14, respectively, 3%, 11%, 24% and 40% 
error^. The fact that a 3% of error in C'^ norm may lead to an order of magnitude error in energy 
(1375%) requires an explanation. 

The reason is the inability of the norm to apprehend "all" the aspects of the similarity 
between two functions. A Sobolev norm [TAY 97], which is also sensitive to derivatives (of arbitrary 
order) may be more appropriate. As a matter of fact, the energy of the free particle in a box may 
be considered to be a Sobolev-type norm: 



n 9 

The mentioned error is obtained multiplying by 100 the £ norm of the difference between both functions, 
being normalized both of them. 



4.5. Warmups for the PRC 
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The Kadanoff wave-functions are smooth when considered in their own "blocked" space, but 
they are not when extrapolated to a refined space. These ideas may lead to further approximations, 
as it is shown in the next chapter. 

The Kadanoff blocking warmup technique is directly generalizable to two or three dimensions. 
Figure 15 shows the wave-functions which result from the warmup of the particle in a bidimensional 
box with fixed b.c. 




Figure 15. 2D Kadanoff warmup for a free particle in a box with fixed b.c. The reduced Hilbert space has 16 
degrees of freedom (4x4 blocks). 

Albeit the Kadanoff warmup has been employed in all the cases exposed in the precedent 
sections, it is not the only possibility. 

Warmups inspired in Wavelets. 

Wavelet theory is a rather developed branch of applied mathematics [LEM 89] [PTVF 97]. We 
have found in it inspiration for the development of a warmup suitable for PRC Similar computa- 
tions have been carried out by S.R. White himself in his "orthlets" theory [WHI99]. 

The general idea is to expand into a multi-resolution or multi-scale basis. This basis contains 
also N states which are organized in a hierarchical way. Figure 16 shows graphically a ID example 
with N = 8. 

The numerical procedure consists of the following steps: 

• The original lattice size is chosen, which must be a power of 2. The technique is, nonetheless, 
generalizable to other values. 

• Establishment of the functional basis. We choose a representation such as that of figure 16. 
The set of functions W(l,i) is spanned by two indices: I means the "level" and i the position. 
Function W(0,0) is uniform. All the rest have zero "average" and support of magnitude 
1\| . 2~^~^^. Functions at the same level have fully disjoint supports. Notice that the set 
{W(l,i)} with I G [0, . . . ,log2(N)] makes up an orthogonal basis of the full original Hilbert 
space. 

• Initialization. We establish ra initial functions: the first ra functions of lowest order in the set 
{W(l,i)}. 

• Iteration. At each step a "superblock" is formed with the m functions which make up the 
approximation so far and a new one taken from the set {W(l,i)}. The superblock hamiltonian 
is diagonalized and we retain just the first ra states. These states make up the approximation 
for the next step. 
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W{0,0) 



W{1,0) 



W{2,1) 



m I a 

W{3,0) W{3,1) W(3,2) W{3,3) 

Figure 16. Multi-resolution basis inspired on wavelets for N — 8 sites. It might also serve as the starting point 
for a larger hierarchy. At each frame, the abscissa represents the spatial position (ID), and the ordinate stands 
for the wave-function. 




The optimum strategy implies to start at the lowest level available and changing from I to 
I + 1 only when convergence at level I has been reached. 

The practical implementation of the method is slow because of the full computation of the 
matrix elements. Even though, the results encourage for further insight and the search for better 
techniques. A couple of sweeps at each level is enough to assure convincing results. For example, 
for the free particle in a box with N = 256 sites and fixed b.c, 3 sweeps yield the energies shown 
in table 5. 





Ground state 


exc. 


2"'' exc. 


3^^ exc. 


Warmup 


0.000155129 


0.000616833 


0.0014313 


0.00253865 


Exact 


0.000149427 


0.000597684 


0.00134471 


0.00239038 



Table 5. Energies obtained with the wavelet-inspired warmup for the free particle in a box with fixed b.c. and 
N = 256 sites. 



4.6 Blocks Algebra. 

The renewal of RSRG techniques is due in great extent to a much freer usage of the block 
concept, even though K.G. Wilson in his study of the Kondo effect had already employed the site 
by site growth pattern. This section is much more abstract than the previous ones, and intends to 
provide a general view of the subject. 

The methods employed in the whole of this thesis may be simply exposed is terms of a blocks 
algebra, which we now develop. 

Let us consider a system with an underlying graph structure Q. Let B{Q) be the set of the 
sub-graphs or blocks of Q. Let B]^ and be elements of B{Q). The links between sites of B^ and 



4.6 Blocks Algebra. 
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Bi^ which do not belong to any of the two blocks shall be contained by definition by the set L\^i. 
Formally, we may write 

Bk + = B„ U U L^i 

The condition for a set of blocks "P to be a blocks algebra is to be closed under the operation 
+: if B]^ and B; belong to the set, then also does Bi^ + Although it is not strictly required, 
we shall only take non-overlapping blocks: addition shall only take place if B]^ n B| = 0. 

A matricial representation of a blocks algebra is obtained by assigning a "block matrix" to 
each B £ V and a "hook matrix" to each pair of blocks L £ V x V. Abusing notation slightly we 
shall denote by Bi< both the block and its associated matrix, and wc shall do the same with hooks. 

The dimensions of the matrices must fulfill the following requirement. A block matrix Bi^ 
is a square matrix of arbitrary dimension d(k), meanwhile L]cv must be a rectangular matrix of 
dimension d(k) x d(l). Even more, it shall be necessary that, if Bg = Bi^k + B|, d(s) < d(k) + d(l) 
(an analogue of the triangle inequality). 

A set {Bi, . . . , BrL, Li2, Li3, . . . , LrL-i,n} shall be termed a representation of a blocks algebra 
of a quantum-mechanical system if there is an orthonormal set of system states {jilJi)} (where 
k G [1, . . . ,n] denotes the block index and i G [1, . . . d(k)] denotes the state within that block) such 
that 

(BOii = {^\>l\ H (W),^ = <tl;!'| H 

We shall now define the operation + acting on the block matrices of a system representation. 
If, in terms of blocks, B]^ + B; = Bg, then operation + acting on the matrices B]^ and Bj, consists 
of: 

• The building of an effective superblock matrix 

which is a square matrix (d(k) + d(l)) x (d(k) + d(l)). 

• The diagonalization of the aforementioned matrix and the retention of the d(s) lowest energy 
eigenstates (which is always possible due to the triangle restriction d(s) < d(k) + d(l)). 

• Truncation so as to retain an effective matrix for the lowest states. If T is the matrix whose 
columns are the weights of the new global states on the old blocks states (first those of B^^ and 
later those of B;) then 

Bs = Bt, + Bv = T(B„ e Bx)Tt 

Now the non-overlapping requirement may be understood: it guarantees the orthonormality 
of the full set of test states: {|tl)]'>}fii^ ^ {|tl)l>}fiV 

• We must still "adapt" the hooks. We may have to perform this process in an explicit way (i.e.: 
having recourse to the real hamiltonian matrix elements and the total wave-functions) . 



Notice that the operation U is only set-oriented: B^. U means the union of both blocks, but without 
restoring the links between them, which is the key notion. 
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A representation of the blocks algebra may always be built if the set minimum blocks (i.e.: the 
sites) are known. In this case, a series of additions among them provide us with a full representa- 
tions. 

Two blocks may be added in presence of a third block which does not overlap with them. We 
shall consider blocks Bi^ and B^ along with an extra block named Bp. Then, the addition 

Bk +p Bt = Tp(Bk©B^©Bp)Tt 
where +p means "addition in presence of Bp", the direct ternary sum means 

(Bk 1-1^,1 
4 I 
■-k.p "-l.p 

and Tp and Tp mean the first d(Bi^ + B;) eigenstates of the effective superblock matrix with the 
weights corresponding to the states \^^) removed (and consequently re-orthonormalized). 

The most useful case of "sum in presence" is that of sum in total presence: B]^ B;, which 
means that the sum is extended in such a way that Bi^ + B^+^^g^ Bp = Q, where R is any partition 
into blocks of the rest of the system. Thus, the superblock hamiltonian represents the complete 
system. Of course, the symbol +j is not uniquely defined: any partition of R = ^ — Bi^ — B]^ may 
do. 




We shall consider the following dynamics on a representation of the system blocks algebra. At 
each step, the block matrix Bg, which corresponds to Bk + B; is substituted for the matrix resulting 
from Bk +7 B; for some valid representation of +j. 

Let Xq be the original representation of the system, obtained by any kind of warmup, and let 
Xt be the representation after t steps. Calling the described procedure an RSRG step and denoting 
the set {Xt} as an RG trajectory, we observe that the described process is a variational finite-size 
RSRG (either implicit or explicit). CBRG, DMRG and PRG are contained as particular instances. 
A fixed point^^ of this RG is that in which the states of each block are the (orthonormalized) 
projections over the corresponding sites of the lowest energy eigenstates of the full hamiltonian. 
In this case, the eigenvalues of the superblock matrix are the lowest energies of the system. 



4.7. Towards a PRG Algorithm for Many Body Problems. 

Having acknowledged the interest of quantum-mechanical problems, it is important to remark 
that the original aim of the effective RSRG techniques was the application to many body problems 
in condensed matter and quantum field theories. The initial development of the BRG was mostly 
performed in the study of such problems (see sections 1.4 and 2.1). ID-DMRG has been employed 
in a large amount of cases and with great success to the analysis of problems of highly correlated 
electrons (Heisenberg model, Hubbard, t — J,... see [HAL 99]). The natural question is: may the 
PRG be generalized so as to work on many particle systems? 



The only fixed point, as a matter of fact, if the blocks algebra contains all the minimum blocks -those of a 
single site. 



4.7. Towards a PRG Algorithm for Many Body Problems. 
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The answer is, at the same time, "yes" and "in practice, not until now" . In principle, being 
an exphcit method, it should suffice to choose a complete basis of the Hilbcrt space, write on it 
the total hamiltonian and start to operate. The problem, obviously, is the combinatory explosion 
in the number of states. For a simple model such as ITF (see section 2.1), the number of states is 
2^. A better alternative is required. 

Let us consider as an example the antiferromagnetic spin 1/2 anisotropic Heisenberg model 
(also known as XXZ model). It is a quantum magnetism model quite well known [MAN gi], which 
is analytically solvable only in ID via the Bethe Ansatz [GRS 96]. The model hamiltonian is given 
by: 



H=Y^ [AS?Sf + i [stS7 + S+Sr) 



(i.i) 

In order to carry out the calculations a full computational platform for the analysis of many 
body problems was developed (see appendix B for computation related questions) . 



Sites Addition and Warmup. 



Let us consider the graph which represents the system and let us establish a building list, i.e.: 
an order for the sites to be added in such a way that the graph is never disconnected and, besides, 
the number of created links at each step is minimum. Let us denote that list by {B^j^j^. 

Let m be the number of states which we intend to obtain from the hamiltonian. Let us take 
the first integer tlq such that 2^° > m and let us build exactly the effective hamiltonian matrix 
for the first tlq sites of the list^^. This matrix is diagonalized and we retain the m lowest states. If 
it is observed that, when "cutting" at the ra-th state we have destroyed a multiplct, the necessary 
states are added in order to complete it. Let vtir be the "real" number of retained states. 

Now the system is prepared for the '^modus ponens" or sites addition process. The effective 
matrices are calculated-^^ for all the operators which may be needed to keep on building, storing 
the matrix elements among the rtir retained states. 

Let us call "active edge" of the system, A, to the set of sites which have not completed their 
links yet. We shall store, therefore, the operators S^, S+ and S" for the sites p of that active edge, 
along with the total hamiltonian matrix Hsb for the system. 

A new site to be added is chosen, which is given by q = Btuq+i. We call G(no) the set of sites 
which already belong to the system. The addition of q is performed in the following way: 

• The stored operators are tensorially right-multiplied by the identity on the new site: 

VpG^, V|^e{z,-,+} S^^S^®l2 
where I2 denotes the identity matrix of dimension 2x2. 

• The Pauli matrices for the new site are tensorially left-multiplied by the identity on the rest 
of the system: 

To be precise, because of the natural division of the Hilbert space into sectors with defined , boxed-matrices 
were employed. Calling "box" C^, to the matrix elements set between states in the sectors {S^} = a and (S^) = b, 
these matrices store explicitly only the boxes which contain non-null elements. 

Split into boxes, of course. 
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V^G{z,-,+} S^ = U,®a^^ 
where denotes, logically, the 2x2 matrices associated to each component. 

• The "links" are calculated between the site p and those of its neighbours which already belong 
to G(no): 

Hnew= E ^S;S^ + ^(S^S-+S + S-) 

peN(q)nG(no) 

• This operator shall be given by a matrix (2mT.) x (2mr), which must be summed to the old 
hamiltonian matrix tensorially right-multiplied by the identity on the new site: 

Hsb = ll-Sb <8) I2 + Hnew 

• Hsb is diagonalized and the m lowest states are retained (again taking care not to break 
multiplets). The operators are renormalized with the truncation operator built from these 
states and the process has completed a RG cycle. 

It is necessary to remark that for the anisotropic Heisenberg model the previously stated 
method, being a mixture of the BRG and the site by site growth techniques (used by Wilson on 
the Kondo problem and in the DMRG), works qualitatively well, and with a reasonable degree of 
precision. 

As an example, table 6 shows some results for square lattices with free boundary conditions of 
small sizes. Conserving only ra = 8 states (thus, tiq = 3), wc obtain less than 5% error. The first 
excited states share the same degree of matching, both in their spin structure and in their energy. 



System 






Exact 


PRG 


Linear N 


= 10, A = 


0.125 


-3.5903 


-3.4995 


Linear N 


= 11, A = 


1.0 


-10.736 


-10.655 


3 X 3, A 


= 0.25 




-4.7443 


-4.5486 


3 X 3, A 


= 0.5 




-6.8797 


-6.3389 



Table 6. Results of the PRG warmup for a many body problem on linear and square lattices, compared to the 
exact diagonalization. The b.c. are always free. 



Absence of Site Substraction Algorithm. 

Once the former cycle has been finished, the problem is to make sweepings, so as to approach 
asymptotically the exact solution. The addition of sites is a well understood variational process. 
The PRG cycle works through the substraction (modus tollens) and posterior addition of the site. 
But, how can we perform this removal? We must start by saying that it has not been carried out 
yet. The reasons for this failure are interesting and shall be displayed in this paragraph. 

Let us consider a family of operators which halve the dimension of the Hilbert space by fully 
removing a site: Ep : Ti^ ^ H^~^. Let us choose now any orthonormal set of states 
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The set |4)|) = Ep |4)|) needs not be orthonormal. It would be if each state belonged to a 
sector with well defined. But in this case, the states of the set |4)-) would mix different sectors 
up. These two complications are the gate through which the difficulties enter. 

Let us suppose that we lift the condition that every state has well defined, since it is only 
technically desirable, but not a required condition for the computation (as it is the orthonormality 
of the states). In this case, we must re-orthonormalize the states. It is not too difficult to do, as 
the reader may observe, since the full chapter 4 has successfully dealt with such troubles. 

But the situation is now rather different. Let us consider the matrix 



Ep Ep |cl)j) 



(where E means that the operator acts on its left) which serves to compute the basis change matrix 
G> which re-orthonormalizes through the fast Gram-Schmidt method. Let us apply this matrix 



on the states so as to obtain a new series 
trivially on site p we have: 



4)- ). For all the operators O which act 



Gj,Gi(4)ic|0|*i) 



while the operators acting on the site must be calculated anew. The removal process has successfully 
ended and we may add it variationally again. 

The problem appears in its whole magnitude when the moment to update the operators comes 
before starting the next RG step. As it is logical, the matrix CT'j must be stored for each value of 
p at every step, and they must be updated. How are we to do it? Let us consider C with q ^ p. 

Its matrix elements on the states 



<^^^) would be 



GUp EqEpG||(|5x) 



The problem is that the scalar products between the states EpEq |ct)i), in which two sites have 
been removed instead of one, can not be "deduced" from the knowledge of the matrices C|. with 
s G S indexing the set of sites. It is required, therefore, to store the matrices of mixed scalar 
products: 



Ep Eq EqEp |4)j) 



for every couple p q (without caring about order) . The reader may now suspect that the game goes 
on. Effectively, in order to renormalize these two-sites operators the knowledge is required of the 
three-points operators... and so on ad nauseam. Due to elementary combinatorics we require 



N + 



+ 



+ 



2^-1 



mxm matrices. I.e.: a quantity of information superior to the storage of the full wave-functions. If 
this difficulty reminds of the closure problem in turbulence theory, the absence of solution reinforces 
the analogy. 



In the end, the reason of the failure of the many body PRG is the same reason for which DMRG 
needed to become explicit when dealing with long-range or multidimensional problems. The open 
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question is: "is it really necessary for an RSRG variational algorithm on non-tree systems to be 
explicit?" 



4.8. Application of the PRG on a Model of Excitons with Disorder. 

Localization, Delocalization and Disorder. 

It is well known since the late 50's [AND 58] [MT61], that the eigenstates of a ID quantum- 
mechanical system with noise strongly tend to be exponentially localized, for any amount of noise. 
Anderson and other authors [AALC 79] proved that, under very general conditions, any decorrelated 
noise repeated the same situation in 2D (see [LR 85] and [RMD 00] for further explanations and a 
series of references) . 

The panorama changed drastically when it was discovered that the absence of correlation 
was crucial. Moreover, the high conductivity of certain materials where disorder is present was 
explained by proving the existence of correlations. 

Recently, F. Dommguez-Adame et al. have proved [RMD 00] the existence of delocalized states 
in ID for systems with decorrelated noise, but with long range order. In their computations they 
studied the model 

'H = J2ei\8,){8,\+J2h)\h){5^\ 

i i,j 

The hopping terms are chosen so as they decay with distance according to a power law: 

h^=J\i-]\-'' 

where J is the coupling between nearest neighbours. The noise appears in the values e-t, which are 
taken from a uniform probability distribution on the interval [—A/2, A/2], where A is called the 
"disorder factor". The relevant parameters are, of course, a and A. 

Although this model was not proposed in order to explain any real system, in [RMD 00] some 
possibilities are considered whose physics might be modelled by it (e.g., planar dipolar systems). 

PRG Application. 

In this case, as in that of the excitons on dendrimers, the numerical experiments involve the 
obtention of the lowest energy spectrum of a system with many degrees of freedom. Although in 
this example there is no underlying graph structure, PRG is probably the best possible tool to carry 
out this computation. The analytical-numerical development of the bidimensional analogue of this 
problem was carried out by F. Dommguez-Adame, J. P. Lemaistre, V.A. Malyshev, M.A. Martm- 
Delgado, A. Rodriguez, J. Rodriguez-Laguna and G. Sierra [DLMRSoi]. 

In this section we shall focus only on the adaptation of the PRG to the mentioned calculations, 
leaving the physical conclusions and the theoretical implications for later work. 

Since the noise is diagonal and the non-diagonal elements follow a scaling law, matrix elements 
of the hamiltonian in the canonical basis are stored by reference. This means that, in the ID case, 
a vector D of N components is kept which contains the diagonal elements D| = et and another 
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vector V, whose k-th component has the value Vi< = Jk""^. In the 2D case only the vector V 
changes, whose k-th component is now Vi^ = Jk""^/^. 

In order to compute the matrix element H|j some steps must be taken: 

• Convert indices i and j into lattice coordinates: Xi, yt, Xj and yj. 

• We compute the square of the distance among them: = (yj — yi)^ + (xj — x-t)^. 

• We obtain the d^-th component of the vector: 'V^2 . 

The rest of the procedure fully coincides with the one discussed in the previous sections. The 
warmup used is the Kadanoff blocking and the computation of the hamiltonian matrix elements 
bears no shortcut due to the absence of graph structure. 

The ID results from [RMD oo] are checked by our own computations. Figure 17 shows, for 
N = 1000 sites-^^ in a ID chain with a = 3/2, three ground states for the system, with A = 1, 
A = 8 and A = 30. 













I 






L 




111 











Figure 17. Ground states of three realizations of the ID system with A'^ = 1000 and a = 3/2, respectively with 
A = 1, A = 8 and A = 30. 



We define the Inverse Participation Ratio (IPR) as 

N 

i=l 

When it is computed on a normalized wave-function, it takes the value 1 only if it is a delta 
function (concentrated on a single point), and it is of order 0(1/N) when the state is delocalized. 
Figure 18 shows the value of the logarithm of the IPR when we vary A for N = 1000 and a = 3/2 
in ID. In the right-hand plot we observe the number of sewings which were needed to reach the 
desired precision of one part in 10-^° for the energies. 

The difficulty of the problem is seen to grow in the transition region, which, as it may be 
readily checked in figure 17, corresponds to wave-functions with the most complex structure. Part 
of the following process, which shall not be considered in the present work, shall be to carry out a 
multifractal analysis of these objects. 

The results of the PRG application for the bidimensional case are now introduced. Figure 19 
shows the dependence of the IPR on the linear size of the lattice for a = 3 and A = 5. Figure 20 
depicts a particular realization of the ground state of the system, obtained with L = 70. 



Our reference value for the computations is A = 10^, but for the sake of a proper visualization, it has been 
reduced by an order of magnitude. 
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Figure 18. The outer plot shows the IPR (in a logarithmic scale) against the parameter A for ID chains with 
a = 3/2 and A'^ = 1000. The inner shows the number of PRG sweeps required to reach the desired convergence 
(one part in 10^°). 




Linear lattice size 



Figure 19. Dependence of the IPR on the length of the 2D lattice. 20 samples were obtained for each of the 
values of L but for the last two ones (10 samples), checking that the statistics were robust for all the cases. The 
straight line has slope —2.05 ± 0.05. 
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Figure 20. Bidimensional wave-function in a 70 x 70 lattice for the ground state of a realization of the system 
with Q = 3 and A = 5. 
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Synopsis. 

5.1. Introduction. 

5.2. RSRG Formalism for Field Evolution Prescriptions. 

5.3. Cells Overlapping Truncators. 

5.4. Applications and Numerical Results. 

5.5. Towards a Physical Criterion for Truncation. 

5.6. Bibliography. 

The present chapter introduces the main lines of an RSRG approach to classical field evolution 
equations which is being developed by the author in collaboration with A. Degenhard-^. Albeit 
some interesting results have already been obtained [DRL oiA] [DRL oiB], the exposition is mainly 
focused on the fundamental ideas, remarking the close connection to the rest of the RG techniques 
in this work. 



5.1. Introduction. 

Field evolution equations, whether deterministic or stochastic, are a fundamental tool in many 
branches of physics, e.g. hydrodynamics, surface growth phenomena, optics or self-gravitating 
media. These equations are usually highly non linear, due to the complexity of the studied phe- 
nomena. 

A constant background topic during the previous chapters was the consideration that physical 
measurements always take place within a certain observational frame, which includes not only the 
reference system and the set of observables to be measured, but also the "grain size" (or ultraviolet 
(UV) cutoff) and the "plate size" (or infrared (IR) cutoff). As relativity theory was born from the 
analysis of the changes of reference frames and quantum mechanics from the study of the influence 

Dept. fur Theoretische Physik -Universitat Bielefeld (Germany) and Physics Dept. -Institute of Cancer Re- 
search (London). 
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of the measurement of different observables on a process, renormalization group theory stems from 
the changes suffered by a physical system when modifying the observational frame. 

Thus, when modelling physical phenomena, partial differential equations (PDE)^ can have no 
real physical meaning. In other words: space can not be assumed to constitute a mathematical con- 
tinuum as far as physical theories are concerned. Observational frames must always have attached 
their IR- and UV-cutoffs. The continuous space-time hypothesis is mathematically equivalent to 
taking the limit in which UV-cutoff 00, which is highly non-trivial. 

As a matter of fact, it may be argued that the continuum hypothesis has been a handicap in the 
conceptual development of physical theories by spawning them with UV and/or IR catastrophes. 
Chaotic dynamics, scaling anomalies, divergences in QFT... were poorly understood pathologies 
before they received the appropriate treatment. 

A field evolution prescription must be an algorithm which provides, once a set of observables 
has been given (along with their uncertainty), the probability density for the same observables after 
some time. Because of overusage of the term "discretization", we have chosen the term aspect to 
denote any such a set of observables. A prescription, therefore, evolves a field aspect in time. 

PDEs may often be converted into field evolution prescriptions, but not always. Models with 
high Lyapounov exponents (for which sometimes not even an existence theorem for the solutions is 
available) yield important problems in practice. For example, currently there is no existence and 
uniqueness theorem for solutions of the 3D Navier-Stokes equation (used to describe newtonian 
fluid mechanics). It might be not a coincidence that there is neither a practical scheme for the 
simulation of high Reynolds number flows. 

The process through which PDE are converted into field evolution prescriptions is usually 
called "discretization". This process reduces the continuous field equations, which have infinite 
degrees of freedom, into a finite-dimensional model. The regions removed from the phase space are 
expected, in some sense, to be irrelevant. In some cases there is a clear-cut algorithm to choose 
the relevant degrees of freedom, but more often (e.g. almost all non-linear PDE) it does not exist, 
and physical intuition is the only guide. 

It is possible to reformulate the question in more general terms if an infinite number of degrees 
of freedom is not assumed for the initial equation. Let us consider the reduction problem from a 
field evolution prescription into another one which has a smaller number of degrees of freedom. 
Let us denote by cj) an aspect associated to the original prescription, and (f)' a renormalized aspect 
associated to the new one. There must be an operator R such that 

= R(*) 

which we shall denote as coarse-graining or truncation operator. Of course, if there is a real 
reduction of the number of degrees of freedom, it shall be impossible to reconstruct the original 
ones, so the inverse operator R~-^ does not exist. 

The evolution prescription for the initial discretization may be given formally by 

9t* = H(l) 

(where H is the evolution generator, for which linearity is not assumed) and the evolution of the 
renormalized aspect shall be given by 



We shall adopt the term PDE with a rather wide meaning, covering any field evolution equation on a 
continuous space-time. 



5.2. RSRG Formalism for Field Evolution Prescriptions. 



121 



at4)' = H'4)' = RHR^4)' 

where is a suitably chosen pseudo-inverse operator. (Of course, RHR^ is a formal operation, 
whose exact implementation depends on the nature of H). The action of different compatible R 
operators constitute a renormalization semi-group. Usual discretizations from continuous space 
into any finite-dimensional aspect are implemented by an R operator with an infinite degree of 
reduction. 

To sum up, we have introduced the following concepts: 

• Aspect of a field: any finite set of observables, preferably along with their uncertainties or 
probability distributions, for an extended object (a field). 

• Field evolution prescription: a set of rules to compute the time evolution of an aspect of a 
field. 

• Truncation or coarse-graining operator: operator which reduces the number of degrees of 
freedom of an aspect of a field. 

• Embedding operator: Pseudoinverse of the former operator, which "refines" an aspect. 

Of course, we might have conserved the more classical terms "discretization" and "discretized 
evolution equation", but we have chosen not to do since the change of meaning might make the 
statements misleading. 

Much interest has been devoted to this subject in the last years, and this work is greatly 
influenced by the recent contributions of N. Goldenfeld and his group: [GHM 98] and [GHMoo]. 
These works were themselves inspired ultimately by the perfect action concept of Hasenfrantz and 
Niedermayer [HN 94], which intended to remove "lattice artifacts" from lattice gauge theories, and 
that of Katz and Wiese [KW 97] about the application of the same concept to fluid mechanics. 

The former works employed a truncation or coarse-graining operator R created by a blocks 
fusion scheme a la Kadanoff. The most important innovation of the present work is the consider- 
ation of a wider set of possible R operators. Sections 5.2, 5.3 and 5.4 deal in some detail with the 
possibility of a geometric approach based on the overlapping of cells [DRL oiA]. Section 5.5 drafts 
the possibility to use a physical criterion (i.e.: problem dependent) in order to choose R [DRL oiB]. 
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Evolution Prescriptions on Partitions. 
Let us consider the following evolution prescription 

9t4>i = H|j(t)j [1] 

This equation may represent any discretization of a linear PDE, whether the algorithm is explicit or 
implicit. The form [1] may also accommodate non-local linear equations and those of higher order 
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in time, along with some complex boundary conditions. Operator H shall be known as evolution 
generator or hamiltonian. 

Some non-linear equations may enter this formalism rather easily. For example, any quadratic 
evolution generator might be implemented as 

9t4^i = Qijk45j45k + Hijct)j [2] 

Surface growing phenomena as governed by the Kardar-Parisi-Zhang (KPZ) equation [KPZ lo] 
[BS 95] or one-dimensional turbulence as described by Burgers equation [BUR 74] [TAF 81] equation 
may be studied this way. 

Truncation Operator. 

The natural interest in the construction of R maps can stem either from limited computer 
resources or due to theoretical reasons, since it may constitute a tool for successively integrating 
out the irrelevant degrees of freedom. 

Field aspects described by equation [1] dwell naturally in a vector space E'^. A truncation 
operator R : E*^ defines a "sub-discretization", and we shall denote the effective field 

as 4)' G E'^. Therefore, the new aspect only has M degrees of freedom. Thus, R must have a 
non-trivial kernel. 

The election of the R operator is the key problem. Ideally, it should depend on the problem at 
hand, i.e. on the field equation and the observables to be measured. In the next section a purely 
geometrical approach is introduced which is independent of the physics of the dynamical system, 
but which uses a quasi-static truncation procedure for a careful selection of the relevant degrees of 
freedom. 

The truncation operator shall be forced to be linear^. This enables us to write its action as 

*i = 

Here we denote the transformed aspect components with primes and capital letter indices. 
Embedding Operator. 

If the R operator had a trivial kernel, an inverse operator might be written, R~^ (the embedding 
operator) and the following equation would be exact 

9tc(), = H,jRrji(f); 

I.e.: the evolution of the exact aspect might be found from the values of the truncated one. 
Therefore, one might integrate a prescription which has only got M degrees of freedom through 
equation 

at(|)i = RiiHijR->'j = H'ij(t>'j [!'] 

Linearity on the truncation operator is a highly non-trivial imposition. It forces the non-mixing nature of 
the evolution generators of different "arity", i.e.: a purely quadratic term can not develop a linear term under 
such a RG construction. This is dubious, e.g., for Navier-Stokes equation, where the nonlinear term is supposed 
to develop a contribution to the viscous (linear) term. 



5.2. RSRG Formalism for Field Evolution Prescriptions. 



123 



where H' is the renormalized evolution generator. Having evolved the reduced aspect, the original 
aspect may be found: 

(t)i(t) = Rr/4)Ut) 
We may express this situation by the commutative diagram: 




[3] 



Equation [1'] requires less CPU time and storage capacity than equation [1] when simulated 
on a computer, but this situation is generally impossible: the truncation operator has a non- 
trivial kernel, so there can be no real inverse. Nevertheless, it is possible to find an "optimum" 
pseudoinverse: an operator RT which fulfills the Moore-Penrose conditions (see [GvL 96]). 



RRPR = R RPRRP = 
(RPR)^ = R^R (RR^)t = RRT 

These equations are fulfilled only if is the singular values decomposition (SVD) pseudoinverse of 
R. is an "extrapolation" operator, which takes a (reduced) discretization from E*^ and returns 
a full E^ one. 

Operators R and play the same role as the truncation operators we have met throughout 
this work, but some extra considerations are needed. 

Operator R as a full matrix retains a big amount of spurious information, inasmuch as the only 
relevant datum is its kernel, i.e.: the degrees of freedom which are neglected. If R^ is the SVD 
pseudo-inverse, then RR^" is the identity on E*^, and R^R is a projector on the relevant degrees of 
freedom subpsace of E^ . This operator shall also be known as the reduction operator. 

The information about these degrees of freedom is stored in the rows of R, which may be read 
as vectors spanning the retained subspace. These vectors need not form an orthonormal set (they 
may even be non-independent). But an orthonormalization operation guarantees that W = R^. 
Whenever we make use of an R operator in this work, we shall mean the equivalence class of 
operators of the same dimension sharing its kernel. Formally, 

R e g£{M^N)/IC 

where QC{M. N) represents the set of matrices taking N degrees of freedom into M. and K. stand 
for the set of internal operations in that space which conserve the kernel. 
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Evolution Prescription Transformations. 

With the pseudoinverse R^, the diagram [3] does not commute in general and its "curvature" 
represents the error within the procedure. The renormalized evolution generator is written as 

H'ij = RKHtjRfj [4] 

where indices are kept for clarity. A renormalized quadratic evolution generator might be written 
like this 

Q'ljK = RiiQijkRj'jR^ [5] 

This expression shall be shorthanded as Q' = RQRP so as to unify notation. Higher degree 
operators are possible, of course, but for the sake of simplicity we shall restrict to the first two-^. 

The schedule for all the simulations which shall be presented in the rest of this work is 

• Present an evolution generator H (at most quadratic) and an initial field ct)(0). 

• Calculate the exact evolution and obtain 45(1). 

• Propose a truncation operator R and obtain its pseudoinverse R^ . 

• Compute the renormalized hamiltonian and the truncated initial field: H' = RHR^ and 

4)'(0) = R4)(0). 

• Simulate the renormalized evolution on 4)'(0) and obtain cf)'(t). 

• Compare (t)(t) and R^4)'(t). 

We distinguish between the "real space error" , which is given by the norm of [(f)(t) — R^4)'(t)] 
(a vector from E'^) and the "renormalized space error", which is the C? norm of [Rct5(t) — 45'("t)] 
(a vector in E*^). These two errors need not coincide. It is impossible for the first one to be zero 
for all the functions in a discrete functional space, but not for the second one. In that case, the 
retained degrees of freedom are exactly evolved despite the information loss. We shall speak of 
perfect action. 



5.3. Overlapping-Cells Truncators. 

In the former section an abstract formalism was introduced for the truncation operators. In 
this section a series of concrete rules shall be proposed to build up R operators by geometric 
principles. 



Trascendental evolution generators, such as that of the Sine-Gordon equation, do not fit in this formalism. 



5.3. Overlapping-Cells Truncators. 
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Let V he a partition of a given manifold Ai (possibly with boundary) into the cells {Ct}[^^. 
Let 45 be a scalar field on this region of space. As it was stated above, an aspect of that field shall 
be defined to be any discrete set of values which attempts to represent the whole knowledge of a 
given observer about the field. The aspect associated to a partition is defined by 



4)1 = / dM,4)(x) 

jm 

where d|a. is any measure on A4. 



Let us choose the interval [0, 1] as our manifold (with boundary) and let Vn denote a regular 
partition of that interval into n equal cells, denoted by = [(i — l)/rL,i/rL]. The truncation 
operator r'^^'^ shall be defined by 



= — -tfjvt; — \°\ 



where \x{-) denotes the measure on A^. In geometrical terms Ct and Cj denote the respective cells 
in the original and the destination partitions. The R matrix elements are the overlap fractions 
between the cells on the destination cell: 

j^jyi^N Overlap between cells Ci and Cj. 
^ Measure of Cj. 

The rationale behind this expression may be described with a physical analogy. Let us consider 
the (pi to denote the density of a gas in each cell of the source partition. The walls between cells 
are impenetrable. Now a new set of new walls is established: the ones corresponding to the new 
partition. The old walls are, after that, removed. The gas molecules redistribute uniformly in each 
new cell. The new densities are the values which constitute the transformed field aspect (see 
figure 1). 
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Figure 1. A region of the plane on which two partitions have been drawn. Black lines delimit the old cells 
{A, B...), meanwhile the red ones correspond to the new ones (1, 2,...). Thus, e.g., there would be no R matrix 
element between cells 1 and C, since they do not overlap. On the other hand, the element Ria must be near 
unity. 
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In more mathematical terms the value of is a linear estimate for 

4)1 = / 4)(x)dx 

conserving total mass: = Si "l^i- 

The resulting r'^^'^ operators shall be called sudden truncation operators. Notice that the 
Kadanoff transformations, in which an integer number of cells get fused into a bigger one, are 
included. 

One of the fundamental differences between usual (integer reduction factor) truncation oper- 
ators and the sudden truncation operators defined in [6] is the possibility of applying truncations 
which remove a single degree of freedom. The situation is shown in figure 2. 

1 



1 

Figure 2. The lower partition has just one degree of freedom less than the one above. 



QuAsiSTATic Transformations. 

The composition of sudden truncation operators takes us to the concept of quasistatic or 
adiabatic truncation operators. These are defined by 

qj^M^N _ |^M^M-|-l]^M-|-l^M-|-2 . . . j^N-l^N 

Of course, qR'^^'^ differs greatly from r'^^'^. The "quasistatic" term is suggested by the ther- 
modynamical analogy which was formerly introduced. Reversibility of a process is related to 
quasistaticity, i.e.: proceeding through very small steps and waiting for relaxation between these. 
In a certain sense, we might expect this transformation to be more reversible and therefore better 
suited to our purposes. 

We formulate analytically a single step sudden transformation, upon which the quasistatic 
operators are based, as 



N - 1 , I 

Iteration of this relation leads to a recursion relation fulfilled by the quasistatic operators 



"^^'^ - M + 1 '^^i^ +mTi'^^i+1-^ 

This recursion allows the calculation of the matrices using no matrix products, which therefore 
renders the approach much more efficient in terms of computer time. 

A question of special relevance is: which are the degrees of freedom retained by this truncation? 
Figure 3 plots five rows of the matrix qR^°^^°. 



5.3. Overlapping-Cells Truncators. 
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80 



Figure 3. Degrees of freedom which are retained by the quasistatic truncation operator proceeding from 80 upto 
20 sites. Rows 1, 5, 10, 15 and 20 are depicted. Notice that these "cells" are now overlapping and have sHghtly 
gaussian nature. 

Each row of matrix qRj^^*^, visualized in figure 3, describes the discretization pattern for 
one of the final 20 cells. Although the functions representing the degrees of freedom are now 
overlapping, they conserve a true geometric nature in real-space. It must be noticed that the 
width of the extreme (leftmost and rightmost) cells is smaller than that of the central one. The 
consequence is a better representation of the boundary conditions. 

The overlapping of cells is not a newcomer in RSRG applications [DEG 99]. In physical terms, 
cells overlapping takes into account inter-cells interaction in an improved way, which is the basic 
global target of RSRG techniques. 

Other Possibilities. 

Even simpler truncation matrices are possible. According to the decimation scheme which is 
common in the literature, a cell out of every N/M = f may be retained, where f takes an integer 
value. This truncation scheme is not well suited for our formalism. The reason is that the chosen 
cells are disconnected. In mathematical terms, the corresponding R matrix would be 



which, along with its SVD pseudoinverse, yield a trivial dynamic. A way to hold linearity, though 
losing the Moore-Penrose conditions, is to use an embedding operator specially designed. 



Special care is required in this scheme to include the boundary conditions appropriately. Of course, 
this is not the standard way to implement this technique. 



Rii — ^i.fl 



fi 





j=f(i-i)+i 



On the other hand, truncation processes based on the Fourier transform with an UV cutoff 
enter the formalism perfectly, but they do not correspond to an RSRG scheme. 
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5.4. Applications and Numerical Results. 

This section discusses some concrete numerical applications, both to linear and non-linear 
examples. 

Heat Equation. 

The heat equation is defined by declaring the laplacian to be the evolution generator. It is 
known that the laplacian operator may be sensibly defined on a great variety of spaces [ROS 97], 
even on discrete ones (see [BOL98] and appendix A). 

Our ID interval shall be [0, 1], split into N cells of width Ax = 1/N. The discrete (or combi- 
natorial) laplacian on the linear graph is given by: 

Uj = -26|j +6|^_j|,i i,je{l...N} 
with fixed boundary conditions: Ln = Lmn =2. The evolution equation is therefore 

The first check of the RSRG scheme shall be to take as initial condition a random increments 
function, generated according to the stochastic differences equation 

4)i+i = (^i+w 

with w a random variable with mean zero and uniformly distributed on an interval of width A. 
Using N = 200, M = 20 (a quite severe reduction of a factor 10) and A = 1/4 we obtain the results 
depicted in figure 4. 




0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 



X X 

Figure 4. A random increments initial condition with 200 cells (left) yields, under exact heat equation evolution 
with K = 1/2 along 500 time-steps with At = 5 ■ 10~^, gives the full line on the right. The triangles mark 
the approximation given by the quasistatic transformation with 20 degrees of freedom. The squares follow the 
decimation approach, and the circles represent the sudden approximation (i.e. usual coarse-graining). 



5.4. Applications and Numerical Results. 
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Real Space Error 


Renormalized Space Error 


Quasistatic 


0.53% 


0.29% 


Sudden 


20% 


19% 


Decimation 


13% 


4.7% 



Table 1 . Errors associated to the functions of figure 4. 

According to the explanations given in section 2, table 1 shows the errors'^ for the curves 
plotted in figure 4. 

Errors are usually smaller in renormalized space. The reason is that in real space two sources 
of error are combined, which we shall know as geometric error and dynamical error. The geometric 
error is due to the transformation R itself. If, e.g., the initial condition belongs to the kernel of R, 
we shall evolve the null function and the error shall be fully geometric. On the other hand, the 
dynamical error is completely given by the "non-closure" of the commutative diagram [3]. Thus, 
a function which is initially orthogonal to the kernel of R may generate immediately components 
on this kernel and be subsequently poorly represented. 

In order to examine the relevant scaling laws [BARg6] a discretization of 5(x— 1/2) is used on 
the 200 cells partition, and is normalized by the condition 

N 

i=l 

Under time evolution this peak initial condition becomes a gaussian function whose width follows 
the law 

w(t) « 

which can be proved to be valid also in the discrete case (see appendix C). 

Using the same parameters as in the previous case, we have performed a quasistatic simulation 
of the problem, showing in figure 5 a log-log plot of the width against time. 

The fit of the data for the quasistatic truncation in figure 5, after a brief transient, fit a straight 
line with slope 0.4990 it 0.0001. The data for the exact evolution of the initial aspect yield the 
same accuracy, but the validity range is complete. The sudden approximation saturates at long 
times and gives a result of 0.66 ± 0.01 for the slope just before that happens. In this case, as it is 
proved in appendix C, decimation gives again the correct result. 



See end of section 5.2 for the specifications on the calculation of errors. 
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1 I ■ — ■ — ■ — I 

1 10 100 1000 

log( Time ) 

Figure 5. Log- log plot of the width of the quasistatic simulation of the delta-function against time. The straight 
line has slope 0.5. 



Low Energy States in Quantum Mechanics. 

The formalism presented here allows us to obtain a very accurate approximation to the low en- 
ergy spectrum of any unidimensional quantum mechanical system, which was the main benchmark 
problem for the previous chapters. 

The transformation H — > RHR^^ may also give an effective transformation for a hamiltonian, 
and its eigenvalues shall be upper bounds to the lowest eigenvalues of the real spectrum^, provided 
that the R and R^ are both orthogonal matrices. 

The diagonalization of the quasistatically truncated laplacian yields very precise values. For 
example, if we choose N = 100 and M = 10, we obtain for the spectrum of — L the values shown 
in table 2. 



Exact 

Quasist. 

Sudden 



0.000967435 
0.000967435 
0.00810141 



0.00386881 
0.00386881 
0.0317493 



0.0087013 

0.00870151 

0.0690279 



0.0154603 
0.0154632 

0.116917 



0.0241391 
0.0247183 
0.171537 



0.0347295 
0.0363134 

0.228463 



Table 2. Low energy spectrum of a particle in a box split into 100 discrete cells, calculated through exact 
diagonalization and two effective variational RSRG techniques. The first is based on a quasistatic transformation 
and the second on a sudden one. 



The bad results for the sudden approximation are a bit misleading. For example, the real space 
error (in norm) for the ground state is only an 11%. But the energy is given by the derivatives 
of the wave-function, so a small difference in norm may yield a big error in energy. 

As a matter of fact, the sudden technique for the case of an integer reduction factor coincides 
with the KadanofF warmup (see section 4.6). When analyzing it, the inability of the norm to 



And also, logically, lower bounds for the highest ones. 



5.4. Applications and Numerical Results. 
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distinguish between two functions was discussed. The energy is in this case a Sobolev-like norm 
[TAY97], which may be much more appropriate^. 

The procedure returns as many eigenvalues as degrees of freedom we have retained, 10 in our 
case. But only the 6 first ones are reasonable. The reason for the great accuracy of the first states 
and the error in the last ones is the precision with which the boundaries are represented. Thus, 
a) the boundary conditions are fairly taken into account, but b) there are not many degrees of 
freedom left for the bulk as we would need in order to have 10 full states correctly represented. 

The method works fine not only for the free particle. Some potentials, such as the harmonic 
oscillator and those studied in the previous chapters, have been analyzed and the results have 
similar accuracy as long as the resulting wave-functions are smooth. 

Kardar-Parisi-Zhang Equation. 

The Kardar-Parisi-Zhang (KPZ) equation is widely used as a model of stochastic and deter- 
ministic surface growth [KPZ 86] [BS 95]. We shall use the following version: 

9tcf) = A|V(|)p + kV^cI) 

which represents a surface in which absorption/desorption phenomena take place, but without 
surface diffusion (which would imply a V^cj) term). 

Although this equation is usually employed in a noisy environment (i.e.: the associated 

Langevin equation is the one which is currently solved), we shall only deal with the determin- 
istic version. The squared gradient term shall be implemented through the quadratic operator^ 

Boundary conditions are imposed which use forward and backward derivatives in the extremes. 

The first test evolves a sinusoidal initial condition 4)(x) = sm(47tx) with x € [0, 1]. A slight 
lack of symmetry of the initial discretization (a "lattice artifact" ) develops a large asymmetry in 
the result. The resolution change is 40 20 and 2000 At = 5 • 10~^ timesteps were taken. Figure 
6 shows the results for A = 2 and k = 1/2. Errors are reported in figure 3. 





Real Space Error 


Renormalized Space Error 


Quasist. 


0.5% 


0.2% 


Sudden 


39% 


38% 


Decim. 


15% 


8% 



Table 3. Errors in the evolution of a sinusoidal initial condition under KPZ, corresponding to figure 7. 



4 2 

Physics for two functions of almost equal L norm may be quite different. One of them may have, e.g., a 
high frequency oscillation or Zitterbewegung superimposed on it. 

^ Warning added in proof: This discretization is known not to conserve the universality class of the continuum 
limit. 
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X X 

Figure 6. A sinusoidal initial condition (left) evolves under KPZ equation (right) with the parameters specified 
in the text. 






Exacto 




Cuasist. □ 




Subito 




\,Decim. 







A different test is obtained using a random increments function, as it was done with the heat 
equation. The chosen parameters are the same as the ones for the previous evolution and the 
results are shown in figure 7, with the errors in table 4. 






Real Space Error 


Renormalized Space Error 


Quasist. 


0.23% 


0.23% 


Sudden 


12% 


11% 


Decim. 


9.4% 


6.5% 



Table 4. Errors corresponding to plots in figure 7. 



Numerical experiments were also carried out on Burgers equation (which is analytically related 
to KPZ) and others, with similar results. 



4.5. Towards a Physical Criterion for Truncation. 



133 



4.5. Towards a Physical Criterion for Truncation. 

The success of the geometricahy based truncation operators described in the previous section 
lead us to think that the reduction of degrees of freedom is feasible. Thus, the application of 
truncation operators which were specific to a given physical model might improve considerably the 
quality of the simulations. This approach shall be termed the physical criterion, as opposed to the 
geometric criterion which was our guide in the previous sections. 

The schedule for the development of such physically oriented truncation matrices is still at their 
first steps. In [DRL oiB] a couple of techniques were proposed which shall be briefly commented in 
this section. 

Short-Time Dynamics. 

A first simple idea is the election of a set of states dependent both on the initial condition 
and its initial "movements". Specifically, if the equation is formally given by equation [1] and 
ct)(t = 0) = 4^0, then the set 

5 = [4)o,H4)o,...,H^-^4)o] 

if duly orthonormalized, constitute the rows of an R matrix (and the columns of its pseudoinverse 
RP) which is highly efficient at short times. 

The reason for the election of the set S should be apparent. For a number of timesteps smaller 
than TL, the state 4'(t) is fully contained in that subspace, which can not be ensured for t > uAt, 
although the approximation is generally still correct for a longer time interval. 

The short-time approach is strongly based on the discrete Taylor series and is analogous to the 
Lanczos diagonalization technique [GvL 96] . This idea has not much utility beyond providing an 
example of a physical criterion, since the long-time dynamics is typically much more interesting. 
Quite often, the approximation becomes unstable beyond its applicability region. 

A curious suitable combination consists in taking degrees of freedom from the quasistatic 
truncation along with those of the short-time approximation. For the heat equation figure 5 shows 
that the fit is perfect after a short transient. This transient is fully removed with the combination 
of both approaches. 

Error Minimization 

An approach which appears to be more promising is the minimization of a suitable magnitude 
which measures the "openness" of the commuting diagram [3]. The usage of such measures is 
not new in the literature. The original formulation of S.R. White's DMRG [WHI 92] used a rather 
similar technique to choose the useful degrees of freedom. In a wide sense, this technique generalizes 
DMRG to non-linear systems. 

So as to make the computations feasible, a single degree of freedom shall be removed at each 
step. We shall employ Dirac's notation from quantum mechanics and call |v) the state to be 
removed. Thus, operators R and R^ shall be any pair fulfilling 



RPR = I~P|^) = 1-14)) (4)| 
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the lack of closure of the diagram [3] for a single time step is given by the operator 

£ = {! + AtH) - RP(1 + AtH)R = P,^) + At (HP,^) + P|^)H - (4)| H |(t>) P,^)) [7] 
Let us consider the Frobenius norm of this operator as the error to be minimized: 

Error = 

When H is a selfadjoint operator, equation [7] is rather easily simplified. Let {|iXi)} be the 
basis of eigenstates of operator H, ordered from lowest to highest eigenvalue: 

N 
i=l 

Let us now decompose the desired state into the basis {|ui)}. 

N 

I*) =X]|J-i|^^i) 
i=l 

Substituting into [7] we obtain the following relation: 

£ij = [l + At(Ei + Ej-(H))]^i^^ 

where (H) = (4)|H|4)). We compute the error and minimize it with respect to variations of the 
with the constraint Yli M-i = 1- Taking previously a first order development in At we get 

Error = ^£?^^1 + 2At (H) 

The minimum value of that expression is obtained when (H) is the lowest eigenvalue of H, so 
14)) = |ui). To sum up, the progressive removal of the eigenstates with lowest eigenvalues is the 
optimum truncation, following a physical criterion. 

The case of the non-sclfadjoint linear evolution generator still admits an analytical development 
which reduces the problem to that of finding the ground state of a difi'erent matrix. Let us consider 
the singular values decomposition of the H operator. 

H = ^Ei K)(vi| 

i 

and let us represent the target vector |4)) in both bases: 

i i 

Both sets of coefficients {xi} and are not independent. If Ci, = (ui|vj), Xi — Yl] CijM-j- 
Expanding expression [7] in the same way as that for the selfadjoint case a parallel expression is 
found to first order in At: 



Error = 1 + 2At (H) 



5.6. Bibliography. 
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The difference appears when computing the minimum value for (H): 

(H) = EiXiM-i = Yl CijHiXiXj 
This shows that the lowest eigenstate of the selfadjoint matrix 

Kij = ^(CtjEt + Cj|EO 

is the required solution. 

The case of selfadjoint operators has a greater physical meaning, since it includes the heat 
equation. In this case the physical recipe just asks to remove the lowest eigenvalue eigenstates of 
the laplacian. These are the states which oscillate most strongly and, in the continuous limit, the 
highest Fourier modes. 

The problem of the computation of the optimum state to be removed forces the minimization of 
a function of degree 2n (where n is the order of the equation) in N variables. Thus, only the linear 
equations offer a panorama suitable for analytical calculations. In the case of quadratic equations, 
the formerly cited work [DRL oiB] presents some results obtained with numerical minimization for 
Burgers and KPZ equations, and the results of the geometrical approach are clearly improved. 
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Conclusions and Future Work. 



• Insight has been gained about the cause of the initial failure of the Blocks Renormalization 
Group (BRG), expressing the conclusions in the development of a rather similar technique 
where the defects were removed: the Correlated Blocks Renormalization Group (CBRG) for 
problems in quantum mechanics both in ID and 2D in the presence of an arbitrary potential. 

• The Density Matrix Renormalization Group (DMRG) has been extended so as to work on 
trees, and it was applied to the analysis of the excitonic spectrum of a family of polymeric 
molecules with fractal character: the dendrimers. 

• The required modifications of the DMRG algorithm in order to be applied to multidimensional 
systems have been studied, arriving at the Punctures Renormalization Group (PRG), which 
uses a single block — more natural than the left-right distinction in > ID. The application to 
the analysis of excitons in disordered systems with long range interaction has been described. 
The difficulties on the path towards a PRG for many body problems were analyzed in some 
detail. 

• Techniques based on real space renormalization group have been employed for the effective 
reduction of degrees of freedom in partial differential equations. Concretely, a theoretical 
framework has been exposed for the sub-discretization process and it has been applied to both 
linear and non-linear equations in ID. 



The Real Space Renormalization Group (RSRG) overflows with theoretical questions which 
require more insight and practical applications which wait for development. The number of pub- 
lications on the field, specially since the DMRG was developed, grows ever larger. In the near 
future, we intend to undertake the following research lines: 

• Application of the PRG, in parallel with analytical RG techniques, to various problems of 
excitons in disordered media. 

• Development of the PRG algorithm for many particles. 

• Development of RSRG techniques for partial differential equations in >1D, specially for fluid 
mechanical problems. 
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A. The Laplacian on a Graph. 

B. Computational Issues. 

C. Asymptotics of the Discrete Heat Equation. 

D. Fast Gram-Schmidt Technique. 

E. Scale Dimension Associated to Energy. 



Appendix A. The Laplacian on a Graph. 

In this appendix some elements of graph theory and their connection to physics are exposed. 
The book of Bollobas [BOL 98] is an excellent introduction to the subject. 

Incidence, Adjacence and Laplacian. 

Let ^ be a simple graph, where V(^) denotes the set of vertices and the set of edges, 
with respective sizes Nv and Ng- The graph structure is given by a neighbourhood structure: 
N : V I— > V{V), which assigns to each vertex the set of neighbouring vertices. The only constraint 
for that structure is symmetry: i G N(j) <^=^ j € N(i). 

Let us assign an arbitrary orientation to the graph by giving a travelling order to each edge. 
The incidence matrix is a Nv x Ng matrix with elements B^j in the set { — 1,0,+1} such that 

{+1 if vertex i is the origin of edge j 
— 1 if vertex i is the end of edge j 
otherwise 

A fundamental matrix related to B is the adjacency matrix, defined to be the Nv x Nv matrix 
whose elements Aii are 



The adjacency matrix is relevant for paths combinatorics on a graph, since the number of 
different ways of traveling from site i to site j in n steps is given by [A^)ij. 

The most important matrix associated to a graph may be the combinatorial laplacian, which 
is defined as the Nv x Nv matrix L = B*B. We define the degree of a vertex dt to be its number 
of neighbours (in physics it is usually called the coordination number), and D to be the diagonal 
matrix Nv x Nv whose (i,i)-th entry is just d(i). Then it is easy to prove the following identity: 




+1 if vertex i is connected through an edge to vertex j 
otherwise 
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L = D-A 

where A is the adjacency matrix. Therefore, the elements of the combinatorial laplacian are just 
the number of neighbours in the diagonal and —1 in all (i, j) entries where vertex i is directly 
connected to vertex j. From now on, we shall drop the adjective combinatorial^. 

It is usual to consider Gelds on graphs for the study of the laplacian. We shall consider a field 
on the graph to be any function (|) : V ^ K, where K may be any numerical field (usually real or 
complex numbers). The functional space to which cj) belongs is a vector space of finite dimension, 
so we may oscillate between the terms "field" and "vector" . 

Spectrum of the Adjacency Matrix. 

Let us denote by {m-i} the set of eigenvalues of the adjacency matrix. Let dmax the maximum 
degree of the graph (the highest coordination number) and dmin the minimum degree. We state a 
series of theorems whose proof is not complicated. 

a) The highest eigenvalue fulfills dmin !^ M-max < dmax- 

b) The maximal degree dmax is an eigenvalue iff Q is regular (i.e.: iff all the vertices have the 
same degree). In that case, it is non-degenerate. 

If —dmax is an eigenvalue, then the graph is regular and bipartite. 

d) un is a subgraph of Q , then the whole spectrum of its adjacency matrix lies within [Xmin and 

M-max- 

Spectrum of the Laplacian 

The laplacian is always a selfadjoint positive semidefinite matrix. The proof is straightforward 
from its definition or, as well, from the fact that the quadratic form (4)| L|4>) (usually called the 
lagrangian form) may be written as a sum of squares. 

Let us denote by |5t) the field concentrated on the i-th vertex. The field jcfj) may be expanded 
in that basis as = (^i |6i) and, therefore: 

(cl)|L|cl)) = 5^dic|>?- J2 E 

iev (ij)eE {tj)eE 

where the summation over (i, j) e E denotes all edges of the graph. The following facts are easy 
to prove: 

a) All the eigenvalues of the laplacian are either positive or zero. 

b) If the graph is regular with uniform degree r, the spectrum of the laplacian and that of the 
adjacency matrix are easily related. If A is an eigenvalue of L, then (j. = r — A is an eigenvalue 
of A (and viceversa). 



Notice that the sign is opposed to the usual one, so as to generate a positive semidefinite operator. We hope 
that this notation, usual in graph theory, does not mislead the reader. 



Appendix A. The Laplacian on a Graph. 
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c) The lowest eigenvalue of the laplacian is always zero. The associated eigenvector is the uniform 
field. The reason is that each row of the matrix has as many off-diagonal —1 elements as the 
number in the diagonal marks. Therefore, the field given by (1, 1, 1, . . . , 1) has always zero 
lagrangian. 

d) The second smallest eigenvalue is of outmost importance. Its magnitude refers to the global 
connectivity. We define the vertex connectivity to be the minimum fraction of the vertices 
which it is necessary to remove for the graph to become disconnected. Then, we get the 
following 

Theorem. The vertex connectivity of an incomplete graph Q is never smaller than the second 
lowest eigenvalue of the laplacian. 

We shall not give a complete proof of this theorem, but only an heuristic argument. The 
eigenfunction associated to the second lowest eigenvalue must be orthogonal to the uniform field. 
Therefore, it must have at least two domains of different signs. But, at the same time, it must 
take the minimum possible value for the lagrangian, which grows with nonuniformity (since it is 
the sum of the squares of the field "jumps"). So, the desired field should have just a single "wall", 
which should be as small as possible. The minimum set of vertices which, when removed, separate 
the graph into two parts would be a nice place to locate it. 

Boundary conditions. 

Probably, the most important constraint with physical applications which may be imposed on 
the functional space of fields on a graph is the presence of non-trivial boundary conditions. We 
shall focus on the Dirichlet (fixed) and Neumann (free) types. 

The discrete analogues of boundary conditions is, of course, highly dependent on our concept 
of "boundary". We shall mark a subset of V(^) to be the border set of the graph (the rest being 
called the bulk). 

So as to make things concrete, we shall consider quasi-regular graphs, which consist of many 
vertices with a common degree r but for a few, which have a smaller degree. By adding a certain 
number of extra vertices we may get all the authentic ones to have a homogeneous degree r, while 
all the new vertices have degree 1, making up what we shall denote to be the closure of the graph. 
The set of added vertices shall be known as tack vertices^ and the vertices which are directly 
connected to them shall be the border vertices. 

At this point we may already consider various possibilities for the boundary conditions. 

a) Free boundary conditions. The acceptable fields take the same value at the tack vertices 

and at the border ones. 

Let us consider any of the border vertices, with index i and degree di = r — Vi, which is 
completed up to r through the addition of tack vertices. Denoting by N (i.) the set of bulk neighbours 
and by T(l) the set of neighbouring tack vertices, the condition would be read as Vk G H(l), 4)^ — 
^i. The action of the laplacian on such a field would be 



The word "tack" denotes, among other things, the loose stitches which fasten a piece of cloth to a frame. The 
term "auxiliar sites", which has already been used [DEGgg], has a similar meaning, but some technical details 
prevent us from using it here. 
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(Lfree = Tcjii - cjjj - ^ 4>k = (r - Vi)(l)i - ^ (^j 

jeN(i) iceT(v) jeN(i) 

V ' 

In other words: the contribution of the tack vertices is exactly cancelled out. Therefore, the 
laplacian does not change form. We might say that the meaning of the free boundary conditions is 
the absence of an exterior world linked to the system: we do not need to consider external elements 
to the graph. 

b) Fixed boundary conditions. The acceptable fields take the value zero on the tack vertices. 
Let i be again the index of a border vertex and, with the notation of the former paragraph, the 
laplacian would be 

(Lfixed = r4>i- *i ~ = T^i - Yl *j = 1- - (A \4>))i 

jeN(i) iceT(i) jeN(t) 



thus we find that it may be rewritten as 



Lfixed = rl - A 

The physical meaning of the fixed boundary conditions is that there is an exterior world, but it 
is trivial. Notice also that the laplacian spectrum with fixed b.c. is closely related to that of the 
adjacency matrix. 

c) Mixed boundary conditions. If the value of the field at the tack vertices is neither zero nor 
the same as that at the border vertices, then the laplacian gets the form: 

(1-mixed \^))i = r4)i- Y ~ $Z 
jeN(i) iceH(i) 

V ' 

constant 

The term marked as "constant" is not linear on the set {(IJilieV! so the laplacian receives an 
inhomogeneous term when the tack vertices are removed. 

This more general type of boundary conditions can take into account the immersion of a system 
into a larger one. 



Random walks, Diffusion and Quantum Mechanics on Graphs. 

Random walks on different kinds of spaces constitute a very important part of the basis of 
theoretical physics. In particular, quantum mechanics of spinless particles may be formulated as 
the statistical theory of random walkers in imaginary time [FH 65] [ITZD 89]. 

Let us consider a large number of non-interacting random walkers in a graph. The number of 
particles at each vertex constitutes a field^. Its time evolution is dictated by the master equation: 



To be strict, the interpretation "number of particles" implies positivity. The requirement may be relaxed 
accepting "anti-particles" . 
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oc ^ Particles which enter — ^ Particles which exit oc 4)^ — dt4)t = — (L |4)))t 

ieN(i) 

whenever vertex i belongs to the bulk. Otherwise, then there is a distinction: 

a) Fixed boundary conditions. Tack vertices (walls) absorb any incoming particle. Therefore, 
border vertices have sinks linked to them. This justifies the "zero" values and the relation to 
the adjacency matrix. 

b) Free boundary conditions. Tack vertices (walls) reflect back any incoming particle. Therefore, 
they may be thought not to exist. In mathematical terms, they take the same value as on the 
border vertices, making the gradient across the walls vanish. 

The combinatorial laplacian should not be considered as a mere discrete approximation to the 
continuous operator used for the Schrodinger equation. The formulation should be the other way 
round: the continuum is just an idealization useful for mathematical physics, meanwhile all the real 
data are discrete. Moreover: the structure of space-time might be discrete itself or, more precisely, 
better approximable by discrete structures than by the continous concepts used nowadays. 
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Appendix B. Computational Issues. 

This appendix deals with certain computational aspects of the present thesis work which may 
be usually neglected. 

We wish to remark that, for the reasons exposed in detail at the preface, all the results shown in 
this memory have been obtained with software written by us and with free software. The programs 
were written in C and C++ (compiled with g++ and gcc from GNU) and it was necessary to write 
a series of libraries, which are cited in the rest of the appendix. 

• Linear Algebra Library, (matrix . cc) A C++ library which implements the classes Vector 
and Matrix, in which sum and multiplication operators have a natural meaning. Standard routines 
for the inversion of matrices, Gram-Schmidt orthonormalization, etc. are implemented. 

Especially important is the exact diagonalization routine. It operates in two steps: obtention 
of a tridiagonal matrix similar to the original one by the Householder method and obtention of the 
eigenvalues through a QL algorithm with implicit displacements [PTVF 97] [GvL 96]. 

• Graphs Library, (graph, cc) The generation and easy manipulation of neighbourhood struc- 
tures was performed with a library which implemented the class called Graph with natural opera- 
tions for the access of the vertices, graphical output and computation of adjacency and laplacian 
matrices, among other functions [BOL98]. 

• Operators Description Library, (odl.cc) The operators used for many body problems are 
manipulated as descriptions in a certain language, having the possibility to act on states, adding, 
multiplying and obtaining their matrix elements on any basis of the Hilbert space. 

• Easy X- Window Graphics Library, (easyx . c) The X- Window graphical system is a very 
extended standard both for PCs and workstations [X 96]. Through the usage of this library of easy 
management, based on the Xlib standard, it was possible to make up animations and to obtain 
real time graphical output. 

• PostScript Graphics Library, (easyps . c) On the other hand, PostScript is the standard 
format for graphics printing [ADO 85]. Many of the pictures of the present memory were produced 
by our own programs using this library. 

These libraries, which are by-products of this thesis, shall be put in the public domain under 
the GPL {General Public License) [GNU 01] in a near future, which shall allow its free usage 
and copying. At the same time, some of these works have already yielded some innovations in 
the computational sciences field, as it may be checked at [GDL 01], where the techniques exposed 
in the odl.cc are being applied to the development of a computational language of geometrical 
constructions with pedagogical applications. 
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Appendix C. Asymptotics of the Discrete Heat Equation. 

It is interesting to notice that many classical results of continuous analysis are also exactly 
true as results of discrete analysis. On the other hand, in many cases it is necessary to make some 
changes so as to preserve the analogy (see, e.g., [GKP 89]). 

In chapter 5 a result was used according to which a delta function under the heat equation 
evolves to a gaussian function which fulfills 

(x^) ^ [1] 

Notwithstanding, the calculations made in that chapter referred to the discrete version of the same 
equation, given in ID by-^: 

ct)t(t + At) = <t>i{t) + (4)i-i(t) - 24)t(t) + 4)i+i(t)) [2] 

This appendix proves that the relation [1] is also exact for the finite differences equation [2] under 
certain constraints for the boundary conditions (b.c.) [DRL 01]. 

Let us suppose the function 4^(0) to be normalized according to 

N 
i=l 

and that the b.c. allow for the exact conservation of that magnitude: either N — > 00, or b.c. are 
free or periodical. It is not allowed to fix the b.c. to any set of values. 

The second moment (x^'j is defined in the continuum by 

(^'(t)>c = I d^x2(j)(x,t) 
We shall assume the values 4>q and 4>N+l (the tack vertices, see appendix A) to be well defined and constant. 
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The discrete analogue would be 



N 



i=l 

We now compute the value for the following time-step: 

N 



N 



i=l 



1=1 
KAt 



KAtAx 



(x2(t + At)) = Ax(iAx)24)t(t + At) 

1 

N 

J]i2ct)t_l(t) - 2 J]i24)|(t) + J]i24),+i(t) 



i=l 



1=1 



i=l 



N 



^i2Ax3ct)i(t) 



i=l 



^2(t)> 



The quantity in brackets may be rewritten as 

N 



N 



+ 1)2 - 2i2 + (I - l)2}ct)t(t) = J] 4>i(t) 



i=l 



i=l 



Due to the conditions previously exposed, the quantity m = ^ 4^1 ("t) may not depend on time, so 
we have 

(x2(t + At)) = (x2(t)) + 2KmAt 
Iterating this equation we arrive at 



(x2(t)) = (x2(0)) + 2Kmt 



as we wanted to prove. 
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Appendix D. Fast Gram-Schmidt Technique. 



In many parts of this work it is necessary to re-orthonormalize a set of vectors which were 
orthonormal before some operation on them took place. This operation is in many cases of "local" 
nature, i.e.: it affects only a small fraction of the components of the vectors. Therefore, a full 
Gram-Schmidt technique is an innecessary waste of computational resources. 

Let us suppose the set of vectors to be orthonormal. An operator O has acted on 

them, so we have a new series of states: 

such that the obtention of the dot products on them is a computationally simple task: 

Cij = (ctJilct^j) 

In the typical case, the dot product of two vectors need N multiplications. We shall assume 
that this is not needed, and matrix C^j is given as an input. 

Once matrix Ctj is known, it is possible to obtain a Gram-Schmidt matrix G| so that the 
states 

make up an othonormal set. We describe now the process to obtain such a matrix. 

The Gram-Schmidt (GS) operation consists, briefly, in the removal from the i-th vector of all 
the "contribution" of the previous vectors (indices 1 up to i — 1). This implies that the GS matrix 
has a structure which may be exploited: G| = whenever j > i. Thus, the orthonormal i-th vector 
is computed only from the previous vectors. 

Let us suppose that the orthonormal set { ct^j^jjzj has already been built. Then the new state 

given by 

i*t) = ict)i)- j;(ct^ii4);>|*;> [1] 

is orthogonal to the previous ones, but is not normalized. Using matrix G| we may rewrite expres- 
sion [1] as 



l*t) = l*i)-E E (4)i|Gil4)„)G{|4);) 

j<i k,l<j 

where the second sum only extends up to k, I < j due to the matrix structure of G. Naturally, G-' 
is a number and not an operator, so we may extract the dot product (4'i|4'k) 

m) = 14)1) -EE CikGj^Gi 14);) = 14),) -EE GiGiQi, |4)i) 

j<l k,l<) 1<1 j <i 
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where the last expression was obtained reorganizing the sums. This last expression may be rewritten 
as 



Ki 



with the coefficients g\ given by 



9i={ 1 if I = i [2] 

otherwise 
Let us suppose that the g\ are known for a given value of i, providing therefore the state |4^?). 
In due time, this state shall be normalized to yield the desired state \^[)- The norm of this state 
shall be: 

{c\>m)= Eg5 9lc(*ii*5)= E ^^^^^ 

j,k<ik<i 5,l^<i 

Therefore, the values of G-^ are obtained from those of g| in a very simple way: 

G} = . [3] 

VEj,k<i9}gicCjk 

The calculation procedure for the G] may be now established in an inductive way. The do 
not require previous values of G>: g^^ = 5y The normalization is direct: 

g;- ' 



and the rest G^^ = 0. 

The values for the g^ may always be calculated from those of the G^ with k < I, as it is shown 
in equation [2] . The usage of this equation along with [3] provides the full algorithm, shown in 
pseudocode 1. 



for all i, j: G(i,j)=0 
G(l,l) ^ 1/V(C(1,1)) 
for i:= 2 upto n 
for j : = 1 upto i — 1 
for k, l:= 1 upto j 

G(i,l) ^ G(i,l) - G(j,k)*G(i,l)*C(i,k) 
G(i,l) ^ 1 
norm 

for k, l:= 1 upto i 

norm ^ norm + G(i,k)*G(i,l)*C(k,l) 
for j:= 1 upto i 

G(iJ) G(iJ)/Vnorm 



Pseudocode 1. Algorithm for the obtention of the matrix G'-. 
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At the end of the previous program, matrix G is complete. 

Appendix E. Scale Dimension associated to Energy. 

In section 2.7 the self-replicabihty properties of certain sets of functions were studied, along 
with their relation to CBRG. In this section some results were shown referring to a certain "pseudo- 
fractal dimension" or "scaling dimension for the energy" . It is a definition with full geometric sense, 
but which is not sensu strictu a fractal dimension [MAN 82] [EDG 90]. 

Let {li'i} be the N components of a normalized wave-function. Let T : V*^ h-> V'^/^ be the 
RG (or coarse-graining) transformation which consists of the following steps: 

$i = ^(452i-l + 4>2i) T4) = 

In other words: we consider the vector formed by the local averages of two sites, and after that 
we normalize it. We may generate a family of states in different vector spaces: 

{(t),T4),T2 4)...T-4>} 

where n is log2(N) — 1. In our applications (although it is not required), N shall always be a power 
of 2. 

Now we take the sequence of expected values of the kinetic energy for these functions, always 
with free b.c: 

^ (4)1 (T-)tH 1^14)) 

In the form of an ordinary sum, 

i=2 

I.e.: the squared sum of the derivative. For "smooth" functions, this sequence approximately fulfills 

log(ETa) f« Co 2m 

The reason is the following: for a smooth function, jumps will be increased in a factor 2 when 
averaging (4 when squared). When normalizing a smaller number of components, these shall be 
increased by a further factor \/2 (2 upon squaring, so a factor 8 so far). But the sum extends only 
to half the components, so we must divide by 2 and a factor 4 remains: the scaling exponent is 
1092(4) = 2 as it was announced. 

We shall denote by ^^energetic scaling exponenf'' or ^''pseudo-fractal dimension associated to 
energy'^ the number e which makes the best fit to 



log(Era) ^ Co + evn 
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In practice we obtain, for example, for the first four normalized polynomials {l,x,x^,x^} on 
2048 sites the values of e: x -> 1.90 ± 0.03, x^ 1.87 ± 0.04 and x^ 1.87 ± 0.04^ As it was 
shown in section 2.7, for the fixed point of the eigenstates of the particle in a box with free and 
fixed b.c. the results of table 4 in chapter 2 are obtained. In the first case they are rather near to 
2, meanwhile in the second one they are clearly Ki 0. 

What physical interpretation do these results have? Let us imagine that these wave-functions 
which were obtained as a fixed point were real wave-functions of a particle and let us make the 
inverse path: from big scales to smaller ones. If e ^ 2 or, at least, is clearly incompatible with 
zero, then the energies of the lowest energy states get smaller and smaller when one takes the 
continuum limit. If the exponent were null, that would imply that the energy of the ground state 
might have a finite limit without being localized. This generation of a "gap" might be interpreted 
as a dynamical mass generation. 



Obviously, there is no sense in trying to find it for a constant function, whose energy is always zero. 
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We present in the following pages an index of commented terms or glossary. Some of these 
terms have been introduced in our work, and in that case they are marked with an asterisk. 

Aspect (*) (of a field) Discrete and finite set of observables referred to an extended physical 
system which aspire to represent the whole knowledge of the observer about it. For the — s>quasistatic 
transformation, e.g., these observables are the integral values on the cells of a given partition of 
space (concept more refined than that of discretization). Sections 5.1 to 5.3. 

Basic Scale (*) Scale at which the physical laws explaining a given phenomenon are simple. 
Thus, e.g., ferromagnetism is explained at atomic scale. Macroscopic physics is obtained from it 
by applying — >RG transformations. Section 1.2. 

Blocks Algebra (*) In the meaning used in this work, blocks are just subsets of a graph, 
on which an internal operation known as addition may be defined. It is possible to interpret RSRG 
as a fiow on a matrix representation of that algebra. Section 4.6. 

BRG, Blocks Renormalization Group, Process for the variational computation of the 
lowest energy spectrum of a system using as Ansatz the lowest energy eigenfunctions corresponding 
to smaller blocks. Section 1.4. Ill and chapter 2 (part I). 

CBRG, Correlated Blocks Renormalization Group (*) Modification of the ^BRG 
which takes into account the correlation between blocks through the introduction of the — i-influence 
and —> interaction matrices. Chapter 2 (part II). 

Dendrimers, Polymeric molecules in which from a monomer stem various branches which 
branch at due time all at once, resulting a system of great symmetry. They have lots of applications, 
from materials science to medicine. Its excitonic spectrum is analyzed through — >DMRG applied 
to — i-trees. Chapter 3 (part II). 

DMRG, Density Matrix Renormalization Group, Drastic modification of the ^BRG 
following the idea that the lowest energy states of the blocks need not be the best bricks for building 
the global state. The system is divided into left and right blocks, and the "most probable" states 
for each block are chosen by fitting to a given — >target state, projecting with a density matrix. 
Esp. section 1.4.V and chapter 3. 

Embedding Operator, Operator which projects a renormalized state into a real space one. 
In QM it is forced to be the adjoint of the — ^truncation operator, while in other applications it 
is only required to be its — >SVD pseudoinverse. Since the renormalized state has fewer degrees of 
freedom, the consecutive application of truncation and embedding does not yield the identity, but 
a projector on the relevant degrees of freedom. Esp. sections 2.1 and 5.2. 

Evolution Prescription (*) (for fields) Discrete computational rule for the evolution 
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(either deterministic or stochastic) of an —> aspect of a field. Sections 5.1 and 5.2. 

Implicit RSRG Method (*) Variational -^RSRG algorithm for QM in which the wave- 
functions are not stored, but a given number of matrix elements of chosen operators on them. 
Otherwise, the method is known as explicit. Sections 3.1 and 4.1. 

Influence Matrix (*) In ^CBRG, an additive modification on the block hamiltonian 
provoked by the presence of a neighbouring block. A given block receives one infiuence matrix per 
neighbour. Section 2.3. 

Interaction Matrix (*) In ^CBRG, part of the total hamiltonian of a system which does 
not belong to the hamiltonian of any block, but stays between a given pair of blocks. Section 2.3. 

ITF, IsiNG Model in a Transverse Field, Quantum model for an uniaxial ferromagnet 
in which the z components of the spins tend to align, while a magnetic field in a perpendicular axis 
(x) tend to prevent that alignment. It has a phase transition even in ID and has also an interesting 
analogy with the classical 2D Ising model. It is analyzed with — >BRG in section 2.1. 

Laplacian, In a rather generic way, it is an operator which allows to compute statistically 
the diffusion of a collective of non-interacting particles on a given space. Especially interesting for 
this work is the laplacian on a graph. Appendix A. 

Observational Frame (*) Any physical observation process requires a certain "grain size" 
(lower or UV cutoff) and a "plate size" (higher or IR cutoff). A pack formed by a reference frame 
(Lorentz, Galilei...) along with the "IR cutoff + UV cutoff" and any other required features for 
observation makes up an observational frame. Section 1.2. 

Patch (*) Set of ^punctures. In ^PRG, region of the system which, at a given RG step, is 
fully well represented in the Ansatz. Section 4.2. 

PRG, Punctures Renormalization Group (*) Variational RSRG calculation of the low 
energy spectrum of a quantum mechanical system via the Ansatz of some functions which represent 
a given region of the system (the — >punctures) and others which represent the rest of it. Chapter 
4. 

Puncture (*) In — >PRG, each of the sites which are well determined in the Ansatz at a 
given step, represented by the corresponding delta state |6p). Section 4.2. 

Quadtree, Structure for the addressing of 2D points or the storage of a bidimensional field 
(e.g. an image) through divisions of a square into 2x2 smaller squares in an interative way. Used 
by 2D ^CBRG. Section 2.6. 

QUASISTATIC Transformation (*) Truncation operator which proceeds by cells overlap- 
ping in an iterative way, removing one degree of freedom at each step. Sections 5.3 and 5.4. 

Renormalized Space, Space of a smaller number of degrees of freedom than real space, in 
which each component represents the weight of a state considered to be relevant. Real states may 
travel to renormalized space via the — ^truncation operators. From renormalized space one may 
jump to real space by using the — ^embedding operators. Esp. sections 2.1, 2.5, 5.1 and 5.2. 

RG, Renormalization Group In physical terms, a displacement of the — >observational 
frame along the scales axis. Esp. chapter 1. 

RSRG, Real Space Renormalization Group RG implementation which only uses blocks 
formed according to geometric criteria. The term is used as opposed to the RG methods based on 
Fourier space. Esp. chapter 1. 
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Scaling Exponents, Many physical laws are ruled by scaling laws or power laws, which 
require a situation with a certain degree of invariance under scaling (aka — >RG) transformations. 
The exponents are robust observables, since they tend to be — ^universal. Critical exponents, 
according to which physical magnitudes diverge near a critical point, are a particular case. Esp. 
sections 1.2 and 2.1. 

Self-Replicability (*) , Property of certain sets of functions on an interval (in d dimensions) 
to be correctly approximated using linear combinations of its copies scaled on a factor 1/2*^ and 
situated on the different quadrants of that interval. Useful concept for the analysis of — >CBRG. 
Section 2.7. 

Sewing (*), Path followed by the — >puncture (or — >patch) through the system in the — S'PRG, 
which must traverse all the links. Esp. section 4.4. 

SVD, Singular Values Decomposition Analogue of the diagonalization for rectangular 
matrices. It allows the easy obtention of a "pseudoinverse" of a given matrix which fulfills the 
Moore-Penrose conditions. Section 5.2. 

Sweep, Process in the second part of traditional — >DMRG, for which one of the two blocks 
(left or right) grows at the expense of the other until an extreme is reached; after that, the process 
continues in the other sense. Analogously, a full path of the free site in the — s-trees DMRG and 
^PRG. Esp. section 3.1. 

Target State, In — >DMRG the block states are chosen as those which fit better to a certain 
state which was obtained for the full system, known as target state. Esp. section 3.1 and 3.2. 

Tree, Connected graph in which there is an only continuous path which does not repeat 
sites between any couple of sites. Especially useful for making calculations with — >DMRG, since 
it allows the use of an — ^implicit method. Section 3.3. 

Truncation Operator, Operator which takes an state from an arbitrary space and returns 
another on a space of smaller dimension, known as the prenormalized space, where each component 
represents the "weight" on a series of states of the original space thought (for some reason) to be 
especially relevant. Esp. sections 2.1 and 5.2. 

Warmup, First step for many RSRG techniques, which initializes the suitable operators before 
the — ^sweeping or — >sewing cycles start. The only requirement is that the matrix elements of the 
operators correspond to real states, but an intelligent election boosts the computation a lot. Esp. 
sections 3.1, 3.2, 3.3 and 4.5. 

Universality Class Set of physical systems which share the same scaling exponents. The 
ubiquity of the phenomenon is explained via the — >Renormalization Group (RG). Section 1.2. 



